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Executive  Summary 


Overview 

Turbulence  modeling  is  a  critical  element  of  modern  CFD  codes  for  flow  analysis.  Models  must  be 
fairly  simple  and  computationally  affordable  to  be  useful  in  repetitive  engineering  calculations,  yet 
they  must  be  sufficiently  representative  of  the  physics  to  capture  the  important  flow  features 
controlled  by  turbulence,  such  as  flow  separation,  skin  friction,  and  heat  and  mass  transfer.  There 
is  particular  need  for  improved  models  that  are  robust  and  accurate  in  rotating  flows,  which  are 
important  in  aerospace  propulsion  systems.  Advanced  aerospace  power  and  propulsion  systems 
employing  magnetohydrodynamics,  and  possible  schemes  for  control  of  electrically  conducting 
external  flow  on  aerospace  vehicles,  also  require  turbulence  models  superior  to  those  presently 
available.  We  have  been  developing  a  new,  structure-based  approach  to  turbulence  modeling 
that  seems  well  suited  towards  meeting  these  needs. 

Under  previous  AFOSR  support,  a  Structure-Based  Model  [1,2]  was  developed  and  shown  to 
perform  better  than  standard  Reynolds  Stress  Transport  models  in  many  complex  flows,  such  as 
flows  with  strong  mean  or  frame  rotation.  Although  the  performance  of  the  original  SBM  was 
considered  quite  satisfactory  [3, 4, 5, 6],  it  was  deemed  necessary  to  simplify  the  approach  into  an 
algebraic  two-equation  model  for  engineering  use. 

This  final  report  describes  the  outcome  of  our  work  towards  developing  a  new  Algebraic 
Structure-Based  model  (ASBM)  for  engineering  use  that  has  been  supported  by  AFSOR  during 
the  last  three  years. 

Stated  Objectives 

The  stated  objectives  of  this  work  program  were  the  following: 

•  to  complete  the  development  of  our  two-equation  algebraic  structure-based  turbulence  model 
(ASBM)  for  complex  engineering  flows 

•  to  evaluate  the  ASBM  for  canonical  and  complex  flows  using  a  wide  range  of  experimental 

•  data  and  data  from  direct  numerical  simulations 

•  to  incorporate  the  ASBM  in  a  typical  engineering  CFD  code  and  make  it  available  for  use  in 
other  CFD  codes 

•  to  extend  our  structure-based  modeling  ideas  to  MHD  flows  and  other  flows  where  structure 
information  is  essential  for  accurate  flow  predictions 

Work  Team 

The  work  was  carried  out  by  a  team  consisting  of  the  Principal  Investigator,  Professor  William  C. 
Reynolds  and  Co-Investigator  Dr.  Stavros  C.  Kassinos.  Drs.  Scot  Haire  and  Carlos  Langer 
completed  their  PhD  dissertations  for  work  carried  under  this  program. 

Results  Obtained 

As  outlined  in  the  attached  technical  report  (Stanford  TF-85)  a  new  Algebraic  Structure-Based 
Model  (ASBM)  has  been  developed.  The  ASBM  has  been  evaluated  in  a  wide  range  of  canonical 
and  complex  flows  using  data  from  both  experiments  and  direct  numerical  simulations.  The 
ASBM  has  been  found  to  produce  excellent  results  for  wall-bounded  flows  with  strong  rotation,  a 
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class  of  flows  with  important  engineering  applications.  The  ASBM  has  been  implemented  in  a 
simple  in-house  CFD  code  at  Stanford,  while  at  the  writing  of  this  report  we  have  initiated 
implementation  of  the  model  in  a  major  commercial  CFD  code.  This  implementation  will  enable 
the  continued  validation  of  the  model  for  a  wide  range  of  complex  flows  and  will  facilitate  its 
dissemination  to  the  industry. 

Future  Plans  and  Synergies 

At  the  writing  of  this  final  report,  we  have  initiated  a  number  of  collaborations  that  aim  at  the  wider 
dissemination  and  use  of  the  results  of  this  effort. 

Under  separate  funding,  one  of  us  (Dr.  Stavros  Kassinos)  has  carried  out  a  series  of  large-scale 
direct  numerical  simulations  of  homogeneous  hydrodynamic  and  magnetohydrodynamic 
turbulence  that  is  sheared  in  rotating  frames  [7],  Results  from  these  simulations  [8]  have  shown 
the  importance  of  using  the  large-scale  enstrophy  equation  as  the  second-scale  equation  in 
turbulence  models  used  for  the  prediction  of  flows  with  strong  frame  rotation.  These  results 
generated  interest  in  the  astrophysical  community  in  using  structure-based  modeling  ideas  in 
investigating  a  range  of  phenomena  including  accretion  discs  and  the  transport  in  solar 
convection  zone  [9].  We  have  also  initiated  a  series  of  collaborations  [10,11,12]  that  aim  at  the 
implementation  of  the  structure-based  model  in  codes  used  for  the  simulation  of  MHD  turbulent 
flows  and  we  hope  to  achieve  dissemination  of  these  ideas  in  that  community  as  well.  Finally,  one 
of  us  (Dr.  Stavros  Kassinos)  has  initiated  collaboration  with  Profs.  John  K.  Eaton  (Stanford 
University)  and  Damian  Rouson  (City  College,  New  York)  aiming  at  the  development  of  structure- 
based  models  for  the  prediction  of  aerosol  transport  in  turbulent  atmospheres.  Finally,  Dr. 
Kassinos  has  established  an  ongoing  collaboration  with  Dr.  Datta  Gaitonde  at  WPAFB  aiming  at 
the  application  of  structure-based  modeling  ideas  in  the  prediction  of  flow  and  heat  transfer 
phenomena  in  advanced  aerospace  power  and  propulsion  schemes  employing 
magnetogasdynamics. 

The  remainder  of  this  Final  Report  is  in  the  form  of  a  Technical  Report  prepared  by  the 
Department  of  Mechanical  Engineering  at  Stanford  University. 
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Abstract 


A  primary  goal  of  RANS  based  modeling  is  to  determine  the  Reynolds  stress  tensor  in 
order  to  close  the  turbulence  problem  at  the  mean  velocity  level.  However,  the  Reynolds 
stresses  alone  are  not  sufficient  to  characterize  complex  turbulent  flows  adequately.  More 
independent  information  on  the  turbulence  is  needed,  especially  in  presence  of  rotation.  For 
example,  complementary  information  can  be  given  by  the  one-point  structure  dimensionality 
tensor.  It  gives  a  measure  of  gradients  in  the  turbulence  —  if  the  turbulence  structures  are 
preferentially  aligned  with  one  direction,  then  the  dimensionality  is  smaller  in  that  direction. 

This  work  uses  hypothetical  turbulent  eddies  to  bring  awareness  of  turbulence  structure 
into  the  turbulence  model.  Averaging  over  a  large  ensemble  of  eddies  produces  a  set  of 
one-point  statistics,  representative  of  the  eddy  field,  and  a  set  of  equations  of  state  relating 
the  Reynolds  stresses  and  the  structure  dimensionality  to  the  eddy  statistics.  An  algebraic 
model  for  the  eddy  statistics  is  constructed  in  terms  of  the  local  mean  deformation  and  two 
turbulent  scales;  the  turbulent  kinetic  energy  and  the  large-scale  enstrophy.  The  algebraic 
model  is  further  sensitized  to  the  presence  of  walls  by  a  blocking  scheme,  which  ensures 
proper  asymptotic  behavior  for  the  Reynolds  stresses  in  the  vicinity  of  walls. 

Contrary  to  existing  ad-hoc  definitions  of  a  second  scale  equation,  the  large-scale  en¬ 
strophy  equation  has  a  fundamental  background;  it  is  derived  from  the  large-scale  vorticity 
equation.  Its  terms  represent  large-scale  processes,  and  their  exact  form  provides  valuable 
guidance  when  making  choices  for  their  closure.  The  exact  forms  also  provide  insight  into 
making  these  model  choices  consistent  with  the  asymptotic  behavior  of  the  exact  terms  in 
the  vicinity  of  walls. 

The  algebraic  model  produces  physically  realistic  Reynolds  stresses  and  structure  ten¬ 
sors  for  different  combinations  of  mean  strain  and  mean  rotation,  with  and  without  frame 
rotation.  The  complete  model,  with  evolution  equations  for  the  turbulent  scales  and  al¬ 
gebraic  equations  for  the  turbulence  structure,  was  successfully  implemented.  Although 
results  for  a  channel  in  a  streamwise  rotating  frame  did  overestimate  stabilization,  the  full 
model  was  found  to  produce  adequate  results  for  a  set  of  channel  flows  in  fixed  frames  and 
in  spanwise-rotating  frames  of  reference. 
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Chapter  1 


Introduction 


1.1  Motivation 

Turbulent  flows  are  present  everywhere,  from  naturally  occurring  flows  to  flows  in  and 
around  mechanical  devices.  Turbulence  is  complex  because  of  the  wide  range  of  temporal 
and  spatial  scales  involved.  The  motion  of  a  simple  viscous  fluid  at  all  scales  is  completely 
described  by  the  Navier-Stokes  equations.  However,  simulation  of  all  scales  poses  a  consid¬ 
erable  challenge  for  computational  resources  now  and  for  the  foreseeable  future. 

The  alternative  to  costly  or  even  currently  out-of-reach  simulations  is  to  perform  some 
form  of  averaging  on  the  turbulent  flow,  drastically  reducing  the  number  of  turbulent  scales. 
Historically  this  started  with  Osborne  Reynolds  and  is  the  engineering  tool  for  mathemat¬ 
ical  analysis  of  turbulent  flows.  The  effect  of  the  discarded  information  is  crucial  in  the 
development  of  the  averaged  flow.  It  enters  the  equations  for  the  mean  momentum  through 
a  new  unknown,  the  Reynolds  stress  tensor,  resulting  in  the  familiar  closure  problem. 

Turbulence  models  aim  to  capture  the  effect  of  the  discarded  information,  on  the  mean 
flow,  by  providing  a  representation  for  the  unknown  terms  in  the  averaged  equations  in 
terms  of  known  quantities.  Traditional  turbulent  models  rely  only  on  the  Reynolds  stresses 
for  bringing  turbulence  anisotropy  information  for  closure. 

However  the  Reynolds  stresses  only  characterize  the  componentality  of  the  turbulence, 
ie.,  which  velocity  components  are  more  energetic.  The  turbulent  field  has  much  more 
information  than  that  contained  in  the  Reynolds  stresses.  The  question  then  is:  what  is  the 
relevant  information  which  should  be  extracted  from  the  original  turbulent  field?  Kassinos 
and  Reynolds  (1994)  demonstrate  that  for  rotating  fluid  turbulence  knowledge  about  the 
turbulence  structure  is  also  crucial.  Among  other  properties,  the  dimensionality  of  the 
flow  is  important.  This  carries  information  about  which  directions  are  favored  by  the  more 
energetic  turbulent  eddies:  if  the  turbulent  eddies  are  preferentially  aligned  with  a  given 
direction,  then  the  dimensionality  is  smaller  along  that  direction. 

The  message  repeatedly  delivered  by  Reynolds  and  coworkers  (1991,  1994,  and  2001)  is 
that  models  based  solely  on  the  Reynolds  stresses  are  not  complete.  Key  illustrations  of 
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the  importance  of  structure  to  capture  flow  physics  are  as  follows: 

•  axisymmetric  expansion  anomaly  (Kassinos,  Reynolds,  and  Rogers  2001),  which  is  be¬ 
yond  the  reach  of  an  equilibrium  algebraic  model  and  hence  not  addressed  here.  Homo¬ 
geneous  turbulence  subjected  to  axisymmetric  expansion  develops  a  larger  Reynolds 
stress  anisotropy  under  slow  deformation  than  under  rapid  deformation.  Moreover, 
when  the  mean  straining  is  removed,  the  stress  anisotropy  may  continue  to  increase. 
This  is  counter  to  the  idea  of  “return  to  isotropy”  built  into  current  day  turbulence 
models,  and  can  be  explained  resorting  to  the  structure  models; 

•  damped  Reynolds-stress  oscillations  for  rapidly  rotating,  unstrained,  initially  aniso¬ 
tropic  homogeneous  turbulence  (Mansour,  Shih,  and  Reynolds  1991),  which  is  beyond 
the  reach  of  an  equilibrium  algebraic  model  and  hence  not  addressed  here.  One-point 
turbulence  models  based  solely  on  the  Reynolds  stresses  for  its  closure  are  unable  to 
predict  the  decay  in  the  Reynolds  stress  anisotropy  in  rapidly-rotating,  unstrained 
turbulence; 

•  material  indifference  of  two-dimensional  (2D)  turbulence  to  rotation  (Speziale  1981), 
which  the  proposed  algebraic  model  does  address.  2D  turbulence  rotated  around  its 
axis  of  independence  is  materially  insensitive  to  the  rotation  (2DMIR).  There  are  only 
kinematic  effects,  so  an  observer  rotating  with  the  frame  or  mean  motion  is  completely 
unaware  of  the  rotation.  Matching  the  2D  material  indifference  limit  may  not  be  as 
important  as  an  awareness  of  the  dimensionality  and  its  degree  of  alignment  with 
the  rotation  vector.  Note  that  a  model  needs  to  be  aware  of  the  dimensionality  of 
the  turbulence  in  order  to  capture  the  2DMIR  effect.  Reynolds  stress  models  are 
only  aware  of  the  componentality  of  the  turbulence.  They  may  be  built  to  satisfy 
properties  of  two-component  turbulence,  e.g.  near-wall  turbulence.  However,  they 
will  not,  in  general,  meet  properties  of  two-dimensional  turbulence,  e.g.  the  2DMIR 
condition.  Ristorcelli,  Lumley,  and  Abid  (1995)  proposed  a  pressure-strain  model 
which  is  consistent  with  the  2DMIR  property  under  certain  conditions.  However, 
their  model  does  not  satisfy  the  2DMIR  property  in  general. 

The  development  of  a  structure-based  turbulence  model  is  a  formidable  task.  Haire 
(2003)  summarizes  the  efforts  of  the  structure-based  turbulence  modeling  group  at  Stanford 
over  the  past  few  years.  Recently  the  effort  has  shifted  towards  the  development  of  an 
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algebraic  structure-based  turbulence  model  (ASBM).  Reynolds,  Kassinos,  Langer,  and  Haire 
(2000)  and  Haire  (2003)  built  an  ASBM  based  on  Kassinos  and  Reynolds  (1994).  It  is 
based  on  an  axisymmetric  eddy-axis  model,  where  the  energy-normalized  Reynolds  stress 
tensor  is  based  on  “axisymmetric”  eddies,  computed  as  in  Kassinos  and  Reynolds.  Eddy 
flattening  was  introduced  in  the  algebraic  model  as  a  way  of  considering  more  general 
energy  distributions  around  the  conceptual  eddies  that  constitute  the  basis  of  the  model. 
The  homogeneous  model  developed  was  found  to  produce  realizable  and  physically  realistic 
turbulent  states  for  diverse  normalized  deformation  rates,  with  and  without  frame  rotation. 
Haire  explored  that  ASBM  in  free  shear  flows,  along  with  a  set  of  scale  equations  ( e.g .  k 
and  s),  and  found  it  produced  adequate  results. 

Flows  of  general  engineering  interest  are  likely  to  be  bounded  by  walls.  The  ASBM 
would  not  be  complete  unless  it  incorporates  some  near-wall  physics,  allowing  it  to  produce 
sensible  solutions  near  walls.  Reynolds  et  al  (2000)  did  propose  an  extension  of  the  ASBM 
for  wall-bounded  flows.  However  coupling  the  complete  model  (with  scale  equations)  to  the 
mean  flow  did  not  produce  satisfactory  results  near  walls. 

The  specification  of  an  ASBM  is  not  complete  without  reference  to  a  set  of  evolution 
equations  for  appropriate  turbulence  scales.  Reynolds,  Langer,  and  Kassinos  (2002)  pro¬ 
posed  a  set  of  scale  equations  based  on  the  large-scale  energy-bearing  motion.  The  develop¬ 
ment  of  the  second  scale  equation  (accompanying  the  turbulence  kinetic  energy  equation) 
is  based  on  the  large-scale  vorticity  equation,  and  its  source  terms  are  representative  of 
large-scale  motion.  It  is  expected  that  coupling  this  equations  with  an  ASBM  can  produce 
a  model  with  superior  predictive  capabilities. 


1.2  Objective 

This  work  sets  out  to  develop  and  assess  a  complete  ASBM  for  wall  bounded  flows,  in 
fixed  and  rotating  frames  of  reference,  revisiting  and  extending  earlier  formulations  of  an 
algebraic  structure-based  model  with  accompanying  scale  equations. 
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1.3  Contributions 

The  most  significant  contributions  of  this  work  are 

•  the  development  of  a  revised  version  of  a  (homogeneous)  ASBM  with  improved  phys¬ 
ical  content,  regarding  eddy  flattening  and  the  specification  of  scalar  parameters  rep¬ 
resentative  of  the  hypothetical  eddy  field; 

•  the  modification  of  the  ASBM  and  the  scale  equations  for  compatibility  with  expected 
near  wall  asymptotic  behavior,  ensuring  the  relative  magnitude  of  the  different  terms 
follows  the  expected  near- wall  trend  of  the  exact  terms; 

•  implementation  of  the  complete  model  for  wall-bounded  flows  (channel  flows)  with  and 
without  frame  rotation,  with  mostly  satisfactory  results.  This  marks  the  first  time 
a  complete  ASBM  has  been  applied  to  wall-bounded  flows,  with  or  without  frame 
rotation. 


Chapter  2 


RANS  and  Turbulence  Structure 


This  chapter  provides  some  essential  background  that  is  needed  in  later  chapters.  Following 
a  brief  review  of  the  equations  governing  the  fluid  motion,  and  of  the  Reynolds  Averaged 
Navier  Stokes  (RANS)  equations,  both  in  a  rotating  frame,  the  closure  problem  is  intro¬ 
duced.  It  is  then  argued  that  in  addition  to  the  Reynolds  stresses,  turbulence  structure 
information  is  needed  to  provide  a  more  complete  description  of  a  turbulent  state,  and  thus 
can  be  used  for  creating  improved  turbulence  models. 

This  provides  the  motivation  for  introducing  a  set  of  turbulence  structure  tensors  and 
their  properties  as  delineated  in  Kassinos  and  Reynolds  (1994)  and  Kassinos,  Reynolds, 
and  Rogers  (2001).  The  structure  tensors  introduced  in  this  chapter  are  then  modeled 
in  Chapter  3  through  algebraic  approximations.  The  same  structure  tensors  reappear  in 
Chapter  4,  where  they  are  used  for  closure  of  the  evolution  equations  for  the  turbulence 
kinetic  energy  and  the  large  scale  turbulent  enstrophy. 


2.1  Equations  of  motion 


The  equations  governing  the  motion  of  an  incompressible  viscous  fluid  in  a  Cartesian  ro¬ 
tating  frame  of  reference  are  given  by  conservation  of  mass  and  momentum  as  in  Greenspan 
(1968), 


duj 

dXi 


=  0, 


dUi  dUi 
dt  JtU^dXj 


+  2  eijk&Uk  =  - 


dP* 

dxi 


&Uj 

v  dxkdxk  ’ 


(2.1) 

(2.2) 


where  w,  is  the  velocity  measured  in  the  coordinate  system  rotating  with  constant  angular 
velocity  f ij,  and  Xj,  t,  P*,  p,  and  v  represent  respectively  the  position  vector,  time,  reduced 
pressure,  density,  and  kinematic  viscosity.  The  reduced  pressure  is  given  by 


l 

P  —  ~  +  ti  —  —  (cjjfcfiJ  Xk){€ipqftpXq) , 
P  ^ 
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where  p  is  the  thermodynamic  pressure,  U  is  such  that  a  conservative  body  force  per  unit 
mass  ( e.g .  gravity)  is  given  by  dU/dx i)  and  the  last  term  represents  the  centripetal  accel¬ 
eration. 


2.1.1  Reynolds-Averaged  Navier-Stokes  (RANS)  equations 

Traditional  Reynolds  decomposition  is  used  in  this  chapter.  The  instantaneous  flow  field 
is  decomposed  into  an  ensemble  average  contribution  and  a  fluctuating  contribution 

Ui=Ui  +  u'u  P  =  p  +  p',  (2.4) 


where  the  overbar  represents  the  averaged  part  and  the  prime  represents  the  fluctuating 
part.  Substituting  into  the  equations  of  motion  and  averaging  produces  the  incompressible 
RANS  equations, 


duj 

dxi 


=  0, 


duj 

dt 


dun  or 

—  +  2eijkn jUk  -  —  +  v 


dP* 


dx 


d2Uj 
dxi  '  "  dxkdxk 


(2.5) 

(2.6) 


2.1.2  Closure 

The  last  term  in  Equation  (2.6)  involves  the  Reynolds  stress  tensor,  Rij  —  u[up  and  is 
representative  of  the  transport  brought  about  by  the  large  scales  of  the  turbulence.  It  is 
an  unknown  in  the  equations  and  the  need  for  its  specification  is  the  penalty  incurred  by 
averaging  the  Navier-Stokes  equations:  the  effect  of  all  the  detailed  information  discarded 
upon  averaging  must  be  accounted  for  in  the  RANS  equations  through  Rij. 

In  traditional  one-point  modeling,  closure  is  obtained  using  only  the  mean  flow,  the 
Reynolds  stresses,  and  a  turbulence  scalar,  i.e .,  any  turbulence  correlation  Qij  (such  as  Rij) 
is  modeled  at  most  as  a  function  of  the  velocity  gradient,  frame  rotation  rate,  the  Reynolds 
stresses  themselves,  and  a  turbulent  scale  like  s,  the  dissipation  rate  of  the  turbulent  kinetic 
energy. 

=  (2-7) 

The  information  carried  within  the  independent  objects  in  Equation  (2.7)  does  not  provide 
any  suggestion  about  the  turbulence  at  different  locations  (nonlocal  information).  In  fact, 
only  the  Reynolds  stresses  describe  the  anisotropy  of  the  turbulence  itself.  However,  the 
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Reynolds  stresses  describe  only  the  componentality  of  the  turbulence.  If  u[  =  0,  then 
Rn  —  0,  indicating  2C  (two  component)  turbulence,  which  can  be  ID  (one  dimensional), 
2D  or  3D.  This  information  about  the  dimensionality  of  the  turbulence  is  not  contained 
within  the  Reynolds  stresses  and  may  be  important  for  modeling.  A  glaring  example  is  the 
property  of  material  frame  indifference  of  2D  turbulence  to  frame  rotation  aligned  with  the 
direction  of  independence,  out  of  the  2D  plane  (Speziale  1981).  In  this  case  the  turbulence 
is  such  that  an  observer  rotating  with  the  mean  motion  will  be  unaware  of  the  rotation 
(material  indifference).  Clearly  a  turbulence  model  needs  to  be  aware  of  the  dimensional 
character  of  the  turbulence  to  exhibit  the  same  property.  Everyday  turbulence  is  not  going 
to  reach  the  2D  limit  under  normal  deformations.  Nevertheless,  as  structures  are  elongated 
in  any  particular  direction  by  the  mean  deformation,  a  decreased  effect  of  the  frame  rotation 
upon  structures  aligned  with  it  should  be  expected. 

Reynolds  (1991),  Kassinos  and  Reynolds  (1994),  and  Kassinos  et  al.  (2001)  repeat¬ 
edly  point  the  shortcomings  of  representing  a  turbulent  flow  solely  through  the  Reynolds 
stresses.  To  incorporate  turbulence  structure,  they  introduce  new  one-point  objects  in  terms 
of  gradients  of  the  turbulent  stream-function  vector,  as  defined  in  the  next  section. 


2.2  The  Turbulent  Stream-Function  Vector  and  the  Turbu¬ 
lence  Structure  Tensors 


The  turbulent  stream- function  vector,  is  introduced  to  explore  and  elaborate  con¬ 
cepts  of  turbulence  structure.  It  is  defined  by 


u. 


:  ^ijk 


dx  7 


M  =  o  ^ 

dxk  ’  dxkdxk  *’ 


(2.8) 


The  turbulence  structure  tensors  are  then  defined  as 


,  r.-MM 

y  dxi  dxj  ’  dxk  dxk  ’  y  dxk  dxj  ’ 


(2.9) 


where  Dij ,  and  Cij  are  respectively  the  structure  dimensionality,  structure  circulicity, 

and  inhomogeneity  tensors.  These  carry  information  about  the  large-scale,  energy-bearing, 
structure  of  the  turbulence  not  conveyed  by  the  Reynolds  stress  tensor. 

Note  that  the  free  indices  in  the  dimensionality  tensor  are  connected  to  the  spatial  deriva¬ 
tives.  It  is  associated  with  the  spatial  isotropy  of  the  turbulence  structures.  If  D\\  —  0, 
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then  the  turbulence  is  instantaneously  2D,  every  is  independent  of  the  x\  direction,  and 
the  turbulence  consists  of  very  long  structures  aligned  with  the  x\  direction. 

The  free  indices  in  the  circulicity  tensor  are  connected  to  the  stream  function  vectors. 
Equation  (2.8)  shows  that  the  stream  function  vector  is  in  turn  related  to  the  fluctuating 
vorticity.  The  circulicity  tensor  is  then  representative  of  the  large  scale  structure  of  the 
vorticity/circulation  field.  If  F\\  =0,  then  there  is  no  large  scale  circulation  about  the 
x\  direction.  Alternatively,  if  Fn  — >  1,  most  of  the  large  scale  circulation  is  concentrated 
around  the  x\  axis. 

The  inhomogeneity  tensor  characterizes  the  degree  of  inhomogeneity  of  the  turbulent 
field.  It  vanishes  in  homogeneous  flows,  which  is  evident  by  rewriting  it  as 


=  MM  =  A  (m Ml) 

u  dxk  dxj  dxk  \Jldxj  J  ' 


(2.10) 


Kassinos,  Reynolds,  and  Rogers  (2001)  note  that  it  is  significant  in  the  near-wall  region 
of  wall-bounded  flows.  However,  away  from  the  wall,  say  y+  >  30,  Kassinos  et  al.  note 
that  the  magnitude  of  Cij  becomes  small  compared  to  the  magnitude  of  the  other  structure 
tensors. 

The  structure  tensors  are  related  to  each  other  and  to  the  Reynolds  stress  tensor.  This 
fundamental  constitutive  relation  is  obtained  by  substituting  the  identity  (Jeffreys  1931) 


trst^opq  —  ^ro^sp^tq  d"  ^rp^sq^to  d"  ^rq^so^tp  ^ro^sq^tp  ^rp^so^tq  ^rq^sp&to  (2.11) 


into  the  Reynolds  stress  definition 

_  Qjj1  Qf{hf 

Rij  —  =  ei9teopq^faT  =  (2.12) 

resulting  in  the  exact  relationship, 

Rij  4-  Dij  -\r  Fij  —  ( Cij  +  Cji )  =  q2 j,  (2.13) 

where  q 2  =  R^k  is  twice  the  turbulent  kinetic  energy. 

In  the  algebraic  structure-based  models  developed  so  far  (including  the  current  effort), 
the  structure  tensors  are  modeled  in  a  (locally)  homogeneous  framework,  where  the  inho¬ 
mogeneity  tensor  vanishes.  The  structure  tensors  are  then  built  such  that  they  obey  the 
fundamental  constitutive  relation,  Equation  (2.13),  without  the  Cij  contribution.  However 
in  near- wall  regions  the  assumption  of  vanishing  C%j  is  poor.  To  reconcile  a  non- vanishing 


Chapter  2.  RANS  and  Turbulence  Structure 


9 


inhomogeneity  tensor  with  the  desired  simplicity  of  the  homogeneous  constitutive  relation, 
the  inhomogeneity  tensor  can  be  absorbed  into  modified  dimensionality  and  circulicity  ten¬ 
sors, 

D%  =  Dij-±(Cij  +  Cji),  (2.14) 

F&  =  Fij-\(.Cij  +  Cji).  (2.15) 

Both  D%j  and  F&  are  symmetric  and  they  equal  their  parent  quantities  for  locally  homo¬ 
geneous  turbulence.  Their  traces  equal  twice  the  turbulent  kinetic  energy, 

D%%  =  F&  =  Rkk  =  q2.  (2.16) 


The  drawback  of  this  apparent  simplicity  is  that  the  diagonal  components  of  the  modified 
structure  tensors  may  become  negative.  In  the  current  models,  and  in  the  model  constructed 
in  Chapter  3,  the  diagonal  components  are  built  such  that  they  are  never  negative. 

Substituting  Equations  (2.14)  and  (2.15)  into  Equation  (2.13),  gives  a  homogeneous-like 
constitutive  relation, 


Rij  +  +  Ffj  =  q2hj- 


(2.17) 


Substituting  the  energy-normalized  structure  tensors, 


_  Rij 


Vh  — 


q 2  5  q*  q* 

into  Equation  (2.17),  and  dividing  by  q2  gives  the  normalized  form  most  commonly  used, 


ICC  _  U 

aij  ~  n2  ’ 


F?9 

fee  _  2J3_ 
Jij 


(2.18) 


rij  +  dK  +  fff^Sij.  (2.19) 

Note  that  only  two  of  the  structure  tensors  are  independent.  Given  the  Reynolds  stresses  and 
the  structure  dimensionality  tensor,  the  structure  circulicity  follows  from  Equation  (2.19). 
Kassinos  and  Reynolds  (1994)  point  out  that  traditional  Reynolds  stress  models,  carrying 
only  the  Reynolds  stress  as  model  variable,  effectively  lump  together  structure  dimension¬ 
ality  and  circulicity  information.  Therefore  they  cannot  isolate  the  correct  dimensionality 
state  of  the  turbulence. 

In  the  next  chapter  an  algebraic  model  for  the  “cc”  tensors  is  developed,  as  no  attempt 
is  made  at  modeling  the  inhomogeneity  tensor.  Modeling  information  for  Cij  is  too  little 
and  too  restricted  to  construct  a  general  model  for  it.  It  turns  out  incorporating  the 
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inhomogeneity  tensor  into  the  other  structure  tensors  should  not  be  critical.  The  uncertainty 
in  the  modeled  structure  tensors,  Dij ,  and  is  likely  large  near  the  wall,  so  even  if  exact 
information  for  C{j  were  available,  it  would  not  improve  the  degree  of  accuracy  of  the  other 
tensors.  More  appropriately,  ignoring  any  explicit  model  for  Cij  does  not  increase  the  error 
in  the  near  wall  area,  as  the  overall  error  in  modeling  the  structure  tensors  should  not 
change  significantly. 


Chapter  3 


Algebraic  Structure-Based  Model 

The  complete  Algebraic  Structure-Based  Model  (ASBM)  consists  of  model  evolution  equa¬ 
tions  for  the  turbulent  kinetic  energy  and  for  the  large  scale  enstrophy  of  the  turbulence, 
along  with  algebraic  model  equations  that  evaluate  the  normalized  structure  of  the  turbu¬ 
lence  (including  the  Reynolds  stresses),  given  the  mean  velocity  gradients,  a  turbulent  time 
scale,  and  a  wall  blocking  parameter. 

This  chapter  describes  the  algebraic  model.  The  objective  is  to  provide  an  algebraic 
procedure  to  evaluate  the  normalized  turbulence  structure  (dimensionality,  circulicity  and 
Reynolds  stresses),  given  the  local  state  of  the  mean  deformation  normalized  by  a  suitable 
turbulence  time  scale.  The  model  is  inspired  by  previous  work  by  Kassinos  and  Reynolds 
(1994)  and  Reynolds  and  Kassinos  (1995)  on  structure  based  turbulence  modeling,  and 
by  recent  work  on  algebraic  structure  based  turbulence  modeling  by  Reynolds,  Kassinos, 
Langer,  and  Haire  (2000)  and  Haire  (2003). 

The  essence  of  the  model  is  the  representation  of  a  general  three-dimensional  (3D), 
three-component  (3C)  turbulent  field  by  the  superposition  of  an  ensemble  of  2D-3C  fields. 
Each  of  these  fundamental  2D-3C  fields  can  be  thought  of  as  an  idealized  turbulent  eddy 
whose  axis,  or  orientation,  corresponds  to  the  direction  of  independence  of  the  2D  field. 
Each  of  these  eddies  supports  motion  around  its  axis  as  well  as  along  it  (3C).  Averaging 
over  a  large  ensemble  of  eddies  gives  relations  for  the  structure  tensors  in  terms  of  statistics 
of  the  eddy  field. 

The  turbulence  structure  and  the  Reynolds  stresses  returned  by  the  algebraic  model 
depend  on  the  choice  of  a  specific  reference  frame,  which  may  or  may  not  be  rotating.  The 
model  is  ideally  applied  to  the  frame  where  the  turbulence  structure  is  at  equilibrium,  or 
closest  to  it. 


3.1  Overview 

Sections  3.2  and  3.3  describe  basic  ideas  associated  with  the  description  of  turbulence 
by  a  field  of  conceptual  turbulent  eddies,  and  the  properties  that  are  assigned  to  each  eddy 
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in  order  to  provide  a  quantitative  “measure”  of  the  turbulence. 

Section  3.4  introduces  ensemble  averaging  over  a  field  of  eddies.  This  provides  statistical 
measures  representative  of  the  turbulence  in  terms  of  properties  of  the  eddy  ensemble. 
Section  3.5  then  finds  constitutive  relations  for  the  Reynolds  stress  and  dimensionality 
tensors  in  terms  of  these  eddy  statistics.  At  this  point  conceiving  a  model  for  the  turbulence 
structure  means  conceiving  a  model  for  the  eddy  statistics  brought  up  in  Section  3.4. 

Vector  and  tensor  representations  of  the  eddy  ensemble  are  modeled  in  terms  of  the  mean 
deformation  in  Sections  3.6  through  3.8.  Finally  in  Section  3.9  the  scalar  eddy  structure 
parameters  are  modeled  in  terms  of  the  local  mean  deformation.  The  model  for  the  scalars 
is  found  by  comparing  states  from  the  model  Reynolds  stresses  and  structure  dimensionality 
tensors  with  target  states,  given  different  mean  deformations. 

Section  3.10  then  sensitizes  this  model  representation  to  the  presence  of  impermeable 
walls,  through  the  use  of  a  blocking  approach  based  on  projection  operators,  aiming  to 
retrieve  the  proper  asymptotic  behavior  for  the  model  Reynolds  stresses. 

Section  3.11  offers  a  brief  summary  of  the  model  equations,  and  Section  3.12  explores  the 
behavior  of  the  homogeneous  (unblocked)  model  to  various  normalized  deformation  rates. 

3.2  Eddy  Concepts 

In  this  model,  individual  eddies  are  the  building  blocks  in  the  description  of  a  turbulent 
field.  In  order  to  proceed  with  the  construction  of  the  model,  tools/measures  for  the  char¬ 
acterization  of  these  individual  eddies  will  be  introduced.  Ensemble  averages  are  then  taken 
over  a  field  of  eddies,  giving  quantities  representative  of  the  eddy  field  along  with  consti¬ 
tutive  equations  relating  the  normalized  Reynolds  stresses  and  the  turbulence  structure  in 
terms  of  the  averaged  eddy  ensemble. 

Three  different  families  of  eddies,  characterized  by  their  associated  motion,  are  con¬ 
sidered;  vortical  eddies,  jetal  eddies,  and  helical  eddies.  To  further  characterize  the  eddy 
motion,  the  possibility  of  eddy  flattening  is  also  introduced. 

Vortical  eddies  represent  2D-2C  fields  where  the  motion  of  blobs  of  fluid  that  compose 
the  eddy  is  around  the  axis  of  the  eddy,  with  zero  velocity  along  the  axis  of  the  eddy.  This 
is  typical  of  irrotational  mean  deformations.  In  fact  the  vortical  mode  is  all  that  is  needed 
to  model  turbulence  in  this  case. 

Jetal  eddies  correspond  to  2D-1C  fields  where  the  motion  is  exclusively  along  the  axis 
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of  the  eddy,  with  zero  velocity  around  the  axis.  Rapid  shear  tends  to  generate  this  kind 
of  motion,  with  a  limiting  state  corresponding  to  counter-flowing  jet  streams  (see  Kassinos 
and  Reynolds  1994). 

Helical  eddies  can  be  thought  of  as  the  superposition  of  vortical  eddies  and  jetal  eddies. 
The  tightness  of  the  helix  will  depend  on  the  relative  contributions  from  the  two  modes. 
Generic  turbulent  fields  are  populated  by  helical  eddies. 

Flattening  is  used  to  describe  the  degree  of  asymmetry  in  the  turbulent  kinetic  energy 
distribution  around  an  eddy.  A  round  eddy  has  no  preferential  direction.  A  fully  flattened 
eddy  may  look  like  a  sheet,  with  all  motion  and  kinetic  energy  concentrated  in  a  plane 
perpendicular  to  the  direction  of  flattening. 


3.3  Eddy  Properties 

To  help  in  the  construction  of  the  different  moments  of  the  turbulence,  each  eddy  is 
envisioned  as  a  conglomerate  of  blobs  of  fluid.  This  approach  follows  from  the  Particle 
Representation  Method  (PRM)  of  Kassinos  and  Reynolds  (1994).  All  blobs  share  a  common 
eddy-axis  direction  vector,  a.  Each  has  its  own  velocity  vector,  V  (with  a  component  aligned 
with  a,  and  the  others  perpendicular  to  it),  and  a  gradient  direction,  n,  perpendicular*  to 
both  a  and  V  (which  are  not  necessarily  orthogonal).  The  PRM  introduces  also  the  stream 
function  vector,  s,  perpendicular  to  V  and  n  (forming  a  mutually  orthogonal  vector  triad).  A 
flattening  of  the  eddy  is  also  allowed;  from  a  round/axisymmetric  distribution  of  its  kinetic 
energy  around  the  eddy  axis  to  a  fully-flattened  distribution,  when  the  kinetic  energy  is 
distributed  in  such  a  form  that  the  eddy  becomes  a  sheet. 

Summarizing,  the  following  properties  provide  a  qualitative  description  of  an  individual 
eddy  in  the  turbulence  field: 

•  eddy-axis  direction,  a; 

•  velocity,  V\ 

•  gradient  direction,  n; 

•  stream  function,  s\ 

•  jetal  and  vortical  motions,  V  —  (Vjet)a  4-  (Vvort)n  x  a ; 

•  flattening  direction,  6,  and 
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•  flattening  intensity. 

An  advantage  of  the  eddy-axis  approach  is  its  relative  ease  in  producing  realizable  states 
for  the  Reynolds  stresses  and  the  turbulence  structure,  while  maintaining  the  constitutive 
relation  rij  +  dtj  +  fij  —  ,  a  consequence  of  the  fact  that  V,  n,  and  s  form  a  vector  triad. 

3.4  Eddy  Ensembles 

Individual  conceptual  eddies  constitute  the  building  blocks  of  the  model.  A  turbulence 
field  is  constructed  by  superposing  a  large  ensemble  of  eddies.  Taking  an  average  over  the 
whole  ensemble  produces  quantities  representative  of  the  field.  Angle  brackets  denote  an 
ensemble-averaged  quantity. 

3.4.1  Eddy-Axis  Tensor 

The  eddy-axis  tensor  is  the  energy- weighted  average  direction-cosine  tensor  of  the  eddy 
axes, 

aij  =  ^2  (V'* aiaj) i  (3«1) 

where  q 2  —  (V2)  is  twice  the  turbulent  kinetic  energy.  The  eddy-axis  tensor  is  determined 
by  the  mean  deformation.  Eddies,  as  material  lines,  tend  to  align  with  the  direction  of 
positive  mean  strain  rate.  They  are  rotated  kinematically  by  mean  or  frame  rotation. 

3.4.2  Eddy-Flattening  Tensor 

The  eddy  flattening  tensor  is  an  ensemble  average  of  the  direction  cosines  of  the  flat¬ 
tening  vector, 

bij  =  ^(V2 bibj ).  (3.2) 

3.4.3  Eddy  Structure  Parameters 
Eddy  Jetting  Parameter,  0 

It  is  representative  of  the  amount  of  energy  in  the  jetal  mode,  and  (1  —  0)  represents 
the  amount  of  energy  in  the  vortical  mode.  Under  ir rotational  deformations,  eddies  remain 
purely  vortical,  0  =  0,  and  in  this  case  the  eddy-axis  tensor  coincides  with  the  circulicity 
tensor.  Rotation  produces  jetal  eddies  (Kassinos  and  Reynolds  1994),  and  in  the  limit  of 
infinite  rapid  distortion  0  — ►  1  for  shear  in  a  fixed  frame; 
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Eddy  Helix  Vector,  7 

The  helix  vector  arises  from  the  correlation  between  the  jetal  and  vortical  modes  of 
motion.  Hence  7^  =  0  for  purely  vortical  or  purely  jetal  turbulence.  Typically  jk  is  aligned 
with  the  total  vorticity  vector  (Kassinos  and  Reynolds  1994).  The  helix  vector  is  the  key 
factor  in  setting  the  shear  stress  level  in  turbulent  fields. 

Eddy  Flattening  Parameter,  x 

If  the  motion  is  not  axisymmetric  around  the  eddy  axis,  the  eddy  is  called  flattened. 
Under  rapid  irrotational  deformation  in  a  fixed  frame  eddies  remain  axisymmetric.  Rotation 
tends  to  flatten  the  eddies  in  planes  perpendicular  to  the  rotation  direction. 


3.5  Reynolds  Stress,  Turbulence  Dimensionality,  and  Turbu¬ 
lence  Circulicity  Constitutive  Relations 

The  quantities  above  describe  the  average  properties  of  a  field  of  hypothetical  eddies. 
Ultimately  they  facilitate  introducing  more  physics,  structure  physics,  into  turbulence  mod¬ 
eling.  The  merit  of  these  quantities  comes  at  the  moment  of  evaluation  of  the  turbulence 
structure  and  the  Reynolds  stresses.  This  evaluation  is  done  by  writing  the  Reynolds 
stresses  and  the  turbulence  structure  in  eddy-axis  coordinates,  through  constitutive  rela¬ 
tions  involving  also  the  structure  parameters.  This  work  departs  with  earlier  structure-based 
approaches  in  the  treatment  of  eddy-flattening.  As  a  result  the  constitutive  relations  depart 
from  the  ones  reported  previously  (Reynolds  et  al  2000,  Haire  2003). 

Conditional  averaging  is  introduced  to  facilitate  the  derivation  of  the  constitutive  rela¬ 
tions.  Conditioned  moments  are  formed  by  averaging  over  eddies  that  share  common  eddy- 
axis  and  eddy-flattening  orientations.  For  an  expanded  derivation  consult  Appendix  A. 
Firstly  an  expression  for  the  dimensionality  tensor  is  proposed.  This  is  then  used  in  the 
constitutive  equation  for  the  Reynolds  stresses.  The  most  general  form  for  the  conditioned 
dimensionality  tensor,  given  an  eddy-axis  direction,  a,  and  a  flattening  direction,  b,  is 

(V^7iiTij\(i)  b) /  V2  =  C\ &ij  +  CzCLi&j  +  Czbibj  +  C±{cLibj  +  <Xjbi)  (3-3) 

H-  C5 (tipqCLpbqbj  -j-  tjpq&pbqbi )  -J-  Cq (tipqdpbq&j  €jpq&pbqCli) , 

where  V2  =  (V2|a,  b).  Applying  the  orthogonality  constraint,  djrij  —  0,  the  normalization 
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constraint,  ( V2ni7ii\d ,  b)/  V2  =  1,  and  setting  C3  =  \  and  C5  =  results  in  the  conditioned 
dimensionality  tensor 

{V2nirij\a, b)/  V 2  =  a{1  -  *[1  -  (anbn)2]}Sij  -  |{1  -  x[l  +  {anbnf\}aiaj  ^ 

+  Xkbj  -  ( anbn)x(aibj  +  ajh). 

Writing  the  velocity  vector  in  an  eddy  coordinate  system,  with  a  component  normal  to 
the  eddy  and  a  component  aligned  with  the  eddy  axis  direction  results  in 

Vi  =  V  ( vmam)ai  +  V  (er5t?;rnsat  )eipqavnq.  (3.5) 

Taking  the  second  moment  of  the  velocity  gives  an  expression  for  the  conditioned  stresses 
in  terms  of  the  eddy-axis  vector  and  the  conditioned  dimensionality  tensor, 

R?/  =  (ViVjlaj)  =  <(umam)2U2a*a#,&>  (3.6) 

+  (( VmQ'7n){€rst'Vr>fl'sQ't)V  ( Cipq&p&jTlq  ~b  ^jpqO'pdiTlq^\(l^b) 

“1“  i^TStVr^s^t)  V  6ipq€jkiCLpTlqQ>kTli\(l1b') . 

Using  Identity  (2.11),  and  defining 

0  =  (Vm)^)  (3.7) 

7s Us  —  (er5t^r^t)(^m^m)^s?  (3.8) 


gives 

Rlf  =  (V24>aiaj)  +  {V2(l  -  4>)[%(5ij  -  ciidj)  +  x{J[l  -  (Onbnftdij  (3.9) 

§[1  d"  {O'n^n)  ] CLi &j  bibj  +  {o>nbn){p>ibj  ~f“  Q'jbi)}]} 

d*  {V2% aitipqQ'pQ'j  H~  CjpqQ'pQ'i)^ {1  —  X\X  ~  {anbn)2]}^qs  +  X^q^s  ~  X(anbn)(asbq)}} • 


Following  Kassinos  and  Reynolds  (1994)  define 

q2rij  =  (Wj),  q2d%  =  {V2ninj),  q2a{j  =  (U2^),  q2b{j  =  (V^bj),  (3.10) 
q2<j>  =  (U20),  q2~fk  =  {V2lk),  q2X  =  (V2X)i  (3.11) 

respectively  the  energy- weighted  normalized  Reynolds  stress  tensor,  eddy-axis  tensor,  eddy¬ 
flattening  tensor,  and  the  energy-weighted  structure  jetal  parameter,  structure  helix  vector, 
and  structure  flattening  parameter. 
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The  final  step  in  the  modeling  process  involves  taking  the  ensemble  average  over  all 
eddies  (eddy- axis  orientations  and  flattening  directions,  a  and  b).  The  ensemble-averaged 
normalized  structure-  and  Reynolds-stress  tensors  become 

=  |[1  —  x(l  ~  anm&mn)]^zj  (3.12) 

—  —  X{\  +  &nmbrnn)\Q>ij  4*  X^ij  ~~  X(ainbnj  4~  Ujn&m)? 

rij  ~  (1  “  4>)  \  ($ij  ~  aij)  +  <t>aij  (3.13) 

H“  (1  $)X  ^nmbmnj^ij  ~~  |(1  4~  ^nm^mn)^ij  ~  bij  4~  Qinbnj  4“  Qjnbni] 

+  (— 1 )(tiprapj  4“  ^jprapi){^[^-  ~~~  X( 1  —  ^nmbmn)]^kr  4"  X^kr  X^fcn^nr}* 

The  normalized  (modified)  circulicity  tensor  is  not  an  independent  quantity.  It  can  be 
obtained  through  the  basic  constitutive  relation  (2.19)  as 

f%  =  $ij-rij-dfj.  (3.14) 


3.6  Eddy-Axis  Tensor  Model 


Following  Reynolds  et  al  (2000),  the  eddy-axis  tensor,  o^,  is  computed  on  the  analysis 
frame,  where  the  turbulence  might  be  at  equilibrium  or  very  close  to  it.  This  is  typically 
a  rotating  frame.  The  eddy  axis  is  computed  with  no  reference  to  the  frame  rotation, 
as  it  is  only  kinematically  rotated  by  it  (Kassinos  and  Reynolds  1994,  Haire  2003).  The 
evaluation  is  divided  in  two  parts.  Initially  a  strained  eddy  axis,  a-,  is  evaluated  based 
on  the  irrotational  part  of  the  mean  deformation.  Next  a  rotation  operation  is  applied, 
sensitizing  the  eddy-axis  tensor  to  mean  rotation.  This  procedure  produces  eddy-axis  tensor 
states  that  mimic  the  limiting  states  produced  under  RDT  for  different  combinations  of 
mean  strain  with  on-plane  mean  rotation,  while  guaranteeing  realizability  of  the  eddy-axis 
tensor.  For  details  see  Appendix  B. 

The  strained  afj  is  given  by 


a  _  1  r  (^ikakj  ^jkaki  3 ^mnanm^ij)T 

°y  =  3^+  ao  +  2A/^L^  ’ 


^2S*kpS>kqa° 


(3.15) 


where  S*j  =  Sij  —  Skk^ij/S  is  the  traceless  strain-rate  tensor,  and  <20  =  1.6  is  a  “slow” 
constant,  ao  is  determined  simultaneously  with  other  model  constants  to  ensure  the  model 
satisfies  a  canonical  boundary-layer  state  (see  Equation  3.68).  Equation  3.15  was  initially 
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proposed  by  Reynolds  et  al.  (2000)  (Haire  2003  proposed  a  derivation  based  on  equilibrium 
considerations).  It  gives  realizable  states  tor  the  eddy-axis  tensor  under  irrotational  defor¬ 
mations.  The  coefficient  of  2  in  the  denominator  allows  it  to  produce  the  correct  limiting 
eddy-axis  tensor  for  arbitrary  rapid  irrotational  deformations. 

The  final  expression  for  the  homogeneous  eddy-axis  tensor,  (for  near-wall  regions  see 
Equation  3.72),  is  obtained  by  applying  a  rotation  transformation  to  the  strained  eddy-axis 
tensor,  afjy 

dij  =  HikHjidfci.  (3.16) 

The  rotation  tensor  is  constructed  as 

+  +  (3.17) 

yjtipp  pp 

where  f i*p  =  ClpqClpq.  Hy  is  required  to  satisfy  the  orthonormal  conditions, 

HikHjk  =  Sij,  HkiHkj  =  Sij.  (3.18) 


The  orthonormality  condition  requires 

hi  =  ijlha  -  Af/2.  (3.19) 

hi  is  determined  with  reference  to  RDT  for  combined  plane  strain  and  rotation,  where 

(  r\  n  n  \ 


(  o  r  o  ^ 


9*  =  - 

«  2 


r  o  o 
ooo 


_  2 


0  ft  0 

-ft  o  o 
ooo 


(3.20) 


By  considering  the  kinematic  distortion  of  line  elements  by  this  mean  flow,  Reynolds  et  al 
(2000)  find  that  the  eddy-axis  tensor  evolves  over  time  to  a  fixed  point  where 


an 


r  +  n 


a  22 


r-n 


a  12  = 


vT2  -  ft5 


for  T  —  ft  >  0. 


2T  ’  2r  ’  2r 

This  determines  h%  for  this  flow.  Results  are  extended  to  generic  flows  with 


ha  =  < 


2  -  2^1(1 +  vT^7)  if  r  <  1 


2-2^1(1-0^17^)  if  r  >  1, 


where 


r  = 


O-pqQqrSrp 

^kn^nmamk 


(3.21) 


(3.22) 


(3.23) 
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3.7  Flattening  Tensor  Model 

The  flattening  tensor  bij  is  modeled  in  terms  of  the  mean  rotation  rate  vector  and  the 
frame  rotation  rate  vector, 

t  (tu+cffifxfij+ofi!) 

ij~  vk  +  cbn{wk  +  cbn{y  [  } 

where  and  ftf  are  respectively  the  mean  vorticity  and  the  frame  vorticity  vectors.  To 
capture  the  correct  sign  of  the  secondary  shear  stresses  for  the  case  of  shear  with  stream- 
wise  frame  rotation  (see  section  3.12.4),  C}>  =  -1.0.  This  choice  for  C&  implies  is  also 
dependent  on  the  particular  choice  for  the  analysis  frame.  This  again  reflects  the  fact  that 
algebraic  models  are  built  for  application  in  frames  where  the  turbulence  is  closest  to  an 
equilibrium  condition. 


3.8  Helix  Vector  Model 


The  helix  vector  7^  is  taken  as  aligned  with  the  total  rotation  vector, 


Ik  =7 


(3.25) 


From  2D  RDT  it  is  found  that  the  shear  stress  correlates  well  with  the  normal  stresses  (see 
discussion  around  Equation  3.45),  suggesting  7  can  be  modeled  as 


7  = 


0 


20(1  ~<t>) 

l+x 


(3.26) 


Equation  (3.26)  introduces  /?,  which  is  indicative  of  the  degree  of  correlation  between  fluc¬ 
tuating  velocities  in  different  directions. 


3.9  Structure  Scalars 

Modeling  0,  /?  (see  Equation  3.26),  and  x  is  a  crucial  part  in  the  construction  of  the 
model.  The  equations  for  these  scalars  are  found  by  comparing  target  turbulent  states  (well- 
defined  Reynolds-stress  and  structure-dimensionality  tensor  configurations),  with  the  result 
obtained  with  the  constitutive  equations  of  Section  3.5,  given  a  local  mean  deformation.  A 


set  of  parameter  values  is  chosen  to  mimic  the  isotropic  turbulent  state  (the  eddy  structure  is 
expected  to  consist  of  axisymmetric  (x  =  0)?  vortical  {<j>  =  0)  eddies).  Tentative  functional 
forms  for  the  same  parameters  are  chosen  with  reference  to  RDT,  by  comparing  RDT 
stresses  and  structure  tensors  with  model  results  (Section  3.5  with  eddy-axis  tensor  and 
eddy-flattening  tensor  computed  for  the  given  deformation).  Interpolation  functions  (along 
with  model  constants)  are  then  chosen  to  bridge  these  limiting  states  (isotropy  and  RDT). 
They  are  selected  specially  to  match  a  canonical  state  of  sheared  turbulence,  observed  in 
the  log  region  of  a  boundary  layer. 

3.9.1  Parameterization 


The  structure  scalars  are  defined  in  terms  of  the  mean  parameters  rjm,  representative  of 
the  ratio  of  mean  rotation  to  mean  strain,  rjf ,  representative  of  the  ratio  of  frame  rotation 
to  mean  strain,  and  a2,  a  measure  of  anisotropy.  The  parameters  are  defined  by 


7]m  — 


(3.27) 


Vf=Vm~  sign(X)y-^, 
2  _ 

ft  —  ClpqCLpq , 


(3.28) 

(3.29) 


where  Aj^r2  is  a  measure  of  the  strength  of  the  mean  rotation,  Q^r2  is  a  measure  of  the 
strength  of  the  total  rotation,  and  52r2  is  a  measure  of  the  strength  of  the  mean  strain. 
These  are  defined  by 


O; 


2  _ 


CLij^lik^kj  >  Clj'  —  QijQfaVtfcj)  S  —  d{j  S]$j ,  X  —  CLijQifcSfcj ,  (3.30) 


and  r,  which  represents  a  suitable  time  scale  of  the  turbulence. 


3.9.2  Exact  Homogeneous  2D  RDT  Solutions 

Throughout  the  model  development  there  is  a  strong  effort  to  make  it  consistent  with 
RDT  solutions,  aiming  to  improve  model  dependability  and  realizability  for  a  wide  range 
of  mean  deformations,  as  well  as  to  obtain  guidance  in  the  functional  shape  chosen  for  the 
structure  parameters.  These  RDT  forms  are  later  blended  with  3C-3D  isotropic  solutions, 
in  such  a  way  that  key  intermediate,  equilibrium  states  are  matched. 
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Reynolds  and  coworkers  found  analytical  solutions  for  limiting  states  of  the  turbulence 
structure  at  infinite  values  of  the  total  shear  (or  strain)  applied,  in  RDT  0f  three-component 
(3C),  two-dimensional  (2D)  turbulence  for  homogeneous  flows  in  the  presence  of  frame 
rotation: 

•  Kassinos  and  Reynolds  (to  be  submitted)  investigated  2D  turbulence  (independent  of 

x\)  in  homogeneous  shear  with  varying  degrees  of  frame  rotation,  aligned 

with  the  mean  flow  vorticity,  ^3  =  — dU/dy  =  — T, 

•  Reynolds,  Kassinos,  and  Langer  (unpublished)  investigated  2D  turbulence  (indepen¬ 
dent  of  £1)  in  homogeneous  hyperbolic  flows  (plane  strain  combined  with  on-plane 
mean  rotation,  e.g .  S12  and  fli2,  £12/^12  <  1),  again  with  varying  degrees  of  frame 
vorticity  orthogonal  to  the  mean  deformation  (aligned  with  the  mean  vorticity  when 
that  is  present). 

These  solutions  identify  a  trapezoidal  region  on  a  (r?/,r/m)  plane  for  which  turbulence 
grows.  The  shear  flows  correspond  to  r/m  =  1  and  any  77/,  while  the  hyperbolic  flows  corre¬ 
spond  to  the  region  0  <  rjm  <  0  and  any  77/.  The  ASBM  attempts  to  capture  the  RDT  value 
of  the  structure  scalars  in  this  region,  in  particular  the  value  of  /?,  which  is  fundamental  in 
determining  the  Reynolds  shear  stress  and  thus  the  production  of  turbulent  kinetic  energy. 
(3  should  decrease  to  zero  outside  of  this  turbulence  growth  trapezoidal  region. 

Choice  of  Analysis  Frame 

To  facilitate  direct  comparison  between  the  limiting  RDT  results  and  the  ASBM,  the 
RDT  results  are  reported  in  a  coordinate  system  aligned  with  the  principal  coordinates  of 
the  eddy-axis  tensor,  where  a2  =  1  (for  these  RDT  limits).  Of  the  flows  analyzed  here, 
in  the  “shear  in  a  rotating  frame”  and  the  “plane  strain  in  a  rotating  frame”  cases,  this 
coordinate  system  coincides  with  that  in  which  the  deformation  is  originally  written  (see 
Equations  3.32  and  3.34).  Note  that  the  mean  deformation  is  given  in  a  frame  rotating  at 
a  rate  such  that  the  mean  deformation  itself,  and  the  limiting  RDT  forms  for  the  Reynolds 
stresses  and  turbulence  structure,  are  constant  in  this  frame. 

For  the  “combined  mean  strain  and  mean  rotation  in  a  fixed  frame”  case  (Equation  3.39), 
the  results  (3.40a-f)  are  quoted  again  in  principal  coordinates  of  the  eddy-axis  tensor.  Note 
however  that  these  are  not  the  same  coordinates  chosen  to  represent  the  mean  deformation. 
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As  the  mean  rotation  grows  compared  to  the  mean  strain,  the  principal  coordinates  of  the 
eddy-axis  tensor  turn,  as  can  be  seen  in  Equation  (3.21). 


Shear  in  a  Rotating  Frame 

Kassinos  and  Reynolds  investigated  incompressible  homogeneous  2D  (independent  of 
xi)  turbulence  subjected  to  mean  shear,  rjm  =  1,  in  a  frame  rotating  about  the  £3  axis, 
where 

Gij  =  TSn6j2,  n{j  =  !Q{ezji,  n{  =  Vf6z3,  (3.31) 

are  respectively  the  mean  velocity  gradient  tensor,  the  frame  rotation-rate  tensor,  and  the 
frame  vorticity  vector.  The  frame  rotation  parameter  is  rjf  =  Qf  /T  =  —  n{2/^i2*  In  matrix 
form, 


(3.32) 


(  0  r  0  ^ 

0 

0 

(  0  -nf  o\ 

0*  _! 

y  -  2 

0 

0 

’  °y  ~  2 

0 

0 

u 

l 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

The  RDT  equations  in 

this  case  are 

9ui  _|_  TV  dU'i  -  -TV*  -19P'  -  OC)f  ,/ 
dt  +  T  2 dx[  ~  l6n  pdxi  Qizz' 


Considering  turbulence  independent  of  the  xi-direction  results  in 


M  _  _r„.o  -  oof  ,/  ^3. --LEE. -oof  u' 

dt  ~  r  2  2“12  2’  dt  ~  Pdx 2  2i  21  11 

du'i  ,  du'2  _  Q 


duo 


1  dp' 


dt 


1  dp' 

pdx3 ’ 


<9xi  <9^2 

Fourier  methods  are  used  to  find  the  solution  for  the  Reynolds  stresses  and  dimensionalities. 
Table  3.1  summarizes  the  limiting  structure  and  Reynolds  stress  tensor  components  for 
infinite  total  shear,  Tt,  for  the  most  general  initial  conditions  considered  by  Kassinos  and 
Reynolds. 

Note  that  in  the  range  0  <  77/  <  1, 


yfruTn 


=  1. 


(3.33) 
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7 If  <  0,77/  >  1 

0  <  r)f  <  1 

O 

II 

77/  =  1 

dn 

0 

0 

0 

0 

d%2 

1/2 

0 

1/4 

1/4 

dsz 

1/2 

1 

3/4 

3/4 

dn 

0 

0 

0 

0 

m 

7//-1 

277/ -1 

1-7?/ 

1 

0 

7*22 

ns 

0 

3/4 

7*33 

0 

0 

1/4 

7*12 

0 

-y/ns(i-ns) 

0 

0 

Table  3.1:  Exact  turbulence  structure  results  for  shear  in  a  rotating  frame.  77/  =  0//r 

Plane  Strain  in  a  Rotating  Frame 

In  the  solution  found  by  Reynolds  et  al  for  hyperbolic  flows  (r)m  <  1),  the  Reynolds 
stress  and  structure  dimensionality  tensors  are  given  by  integrals  of  hypergeometric  func¬ 
tions;  there  is  no  simple  algebraic  form  as  in  the  shear  case  analyzed  above.  To  simplify 
the  model  developed  later,  a  polynomial  approximation  for  the  exact  solution  along  the 
line  rjm  =  0  is  determined.  Then  the  parameters  in  the  region  0  <  r}m  <  1  are  found  by 
interpolating  between  the  results  along  the  lines  rjm  =  0  and  rjm  =  1. 

Here  the  analysis  is  restricted  for  incompressible  homogeneous  2D  (independent  of  x\) 
turbulence  subjected  to  mean  plane  strain  in  a  frame  rotating  about  the  X3  axis.  In  matrix 
form,  the  mean  strain  rate,  mean  rotation  rate,  and  mean  frame  rotation  rate  tensors  are 


O 

O 

t-4 

0 

0 

0 

0 

1 

0 

Stj  = 

O 

1 

O 

,  ilij  — 

0  0  0 

.  04-5 

nf  0  0 

(  0  0  0  j 

0 

0 

0 

O 

O 

O 

where  again, 


ni  =  si's* 3. 


(3.35) 
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Considering  turbulence  independent  of  the  aq-direction,  the  RDT  equations  for  this  case 
are 


M_ri2g  =  _ru;_2nU, 


dn <2 
dt 

dt 

du\ 


■r^  =  ru'2. 


nfu' 

Pdx  2  21  11 


rX2^k  =  _I# 

2  &c2  pdxi 


+  ' 


2  _ 


0. 


Fourier  methods  are  again  used  to  find  the  solution  for  the  Reynolds  stresses  and  dimension¬ 
alities.  The  solution  for  different  tensor  components  involves  integrals  of  hypergeometric 
functions.  Tabulated  values  of  the  solution  (rjm  =  0)  are  given  on  Table  3.2.  The  frame 
rotation  parameter  is  now  Tff  =  Q,f  /T  =  fii2/5n,  the  ratio  of  frame  vorticity  to  the  mean 
strain  rate. 

A  good  fit  to  the  discrete  values  from  Table  3.2  is  given  by 

\2\  61 


dsn  =  —  —  — 
"22 


rn  —  0.145 


0.342^|^  +  (1  -  0.342)  'f ) 


'(2r?/)2 

(Vvf)2\9] 

3/4 

l  3/4  J 

(3.36) 

(3.37) 


for  0  <  (2 r]f)2  <  3/4,  where  rtf  =  f]{2/r.  If  (2t?/)2  >  3/4,  then  r33  =  cfo  =  1. 

Note  that  then  for  plane  strain  in  a  rotating  frame,  each  line  on  Table  3.2  gives 

I»12| 


\/rur22 


■  0.95. 


(3.38) 


Combined  Mean  Strain  and  Mean  Rotation  in  a  Fixed  Frame 


To  help  determine  the  interpolation  function  between  r)rn  =  0  and  r)m  =  1,  the  solution 
for  the  2D  turbulence  problem  of  mean  strain  combined  with  mean  rotation,  with  no  frame 
rotation,  is  determined.  In  a  coordinate  system  x,  ,  the  mean  deformation  has  the  following 
matrix  form, 


o 

o 

1 

1 

2 

o 

o 

j  **ij  ~  2 

o 

o 
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' 

0  fi  0 
-Q  0  0 
0  0  0 


’  ^ij 


/  0  0  0  ^ 


0  0  0 
0  0  0 


(3.39) 


Chapter  3.  Algebraic  Structure-Based  Model 


25 


rn 

V22 

^33 

ri2 

^33 

0 

0.00 

0.50 

0.50 

0.00 

0.50 

0.50 

0.01 

0.000020457 

0.49997045 

0.500009087 

0.0030684 

0.50001743 

0.49998258 

0.1 

0.0020465 

0.497034 

0.5009196 

0.030595 

0.501755 

0.498245 

0.2 

0.0081940 

0.48799 

0.503816 

0.060627 

0.507185 

0.492815 

0.3 

0.018453 

0.47239 

0.509156 

0.089429 

0.51683 

0.48317 

0.4 

0.032777 

0.44929 

0.51793 

0.11607 

0.53179 

0.46821 

0.5 

0.050892 

0.41692 

0.53219 

0.13904 

0.55422 

0.44578 

0.6 

0.071728 

0.37173 

0.55654 

0.15541 

0.58857 

0.41143 

0.7 

0.091148 

0.30529 

0.60356 

0.15810 

0.64600 

0.35400 

0.75 

0.095793 

0.25723 

0.64698 

0.14834 

0.69279 

0.30722 

0.79 

0.091532 

0.20464 

0.70383 

0.12897 

0.74872 

0.25128 

0.8 

0.088294 

0.18829 

0.72341 

0.12140 

0.76701 

0.23299 

0.85 

0.040454 

0.067237 

0.89230 

0.048868 

0.913094 

0.086906 

0.86 

0.017884 

0.028178 

0.953938 

0.021009 

0.963211 

0.036789 

0.865 

0.0033232 

0.0050950 

0.9915818 

0.0038486 

0.9933132 

0.0066868 

0.866 

0.000083844 

0.00012784 

0.99978831 

0.0000968 

0.99983204 

0.00016796 

V^4 

0.00 

0.00 

1.00 

0.00 

1.00 

0.00 

Table  3.2:  Integration  results  for  plane  strain  in  a  rotating  frame.  7]f  =  0//T 


The  RDT  equations  in  this  frame  are 


du[ 
dec ■ 

M  +  (r  +  n),2M  +  (r_n)a:iM 


=  o, 

=  -(r + n)u'26u  -  (T  -  n)«i  si2  -  -  2n{yz. 


The  velocities  and  coordinates  are  first  transformed  into  velocities  and  coordinates  in  a 
Cartesian  frame  aligned  with  the  axis  of  independence  of  the  2D  turbulence,  where  an  =  1. 
To  solve  the  equations,  the  coordinate  system  is  further  transformed  into  one  deforming 
with  the  mean  flow  (Rogallo  transformation).  In  this  system  the  convective  terms,  with 
their  explicit  dependence  on  the  coordinate  axes,  vanish.  Fourier  methods  can  then  be  used 
to  find  the  solution.  In  this  case  the  integration  for  the  structure  dimensionality  and  the 
Reynolds  stress  tensors  can  be  carried  out  analytically  in  the  infinite  deformation  limit, 
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Tt  -►  oo,  to  obtain  (in  a  coordinate  system  aligned  with  the  principal  direction  of  the 
eddy-axis  tensor) 


?’22  =  ^33  — »  - 

ni  -+  f  (1  -ln2)C2 


l  +  -(l-ln2)<2 


T-l 


ri2  -*•  -2(1  -  ln2)C 


1  +  |(1  -  ln2)C2 
l  +  |(l-ln2)C2 


1-1 


1-1 


1-1 


d22  -» r-33  -  ^(7T2  -  12 In 4  +  31n24)C2  |l  +  |(1  -  ln2)<2 
^33  r22  +  ^(8  +  7r2  -  161n4  +  31n24)C2  [l  +  ^(1  -ln2)C2 


) 

1-1 


c2  = 


V m 


1  -r?2  ’ 

x  '/m 


(3.40a) 

(3.40b) 

(3.40c) 

(3.40d) 

(3.40e) 

(3.40f) 


where  rjm  =  Sl/T.  Again,  to  facilitate  comparison  with  the  algebraic  model,  the  structure 
and  Reynolds  stress  tensors  are  written  in  principal  coordinates  of  the  eddy-axis  tensor 
(found  with  resort  to  Equation  3.21),  instead  of  the  original  coordinates  in  which  the  mean 
deformation  was  written  (Equation  3.39). 

Note  that  in  combined  mean  strain  and  mean  rotation  (no  frame  rotation), 


I712!  « 0.96.  (3.41) 

V rn  7*22 


3.9.3  ASBM  Solutions 

To  determine  the  structure  parameters  based  on  the  information  gathered  from  the  exact 
RDT  solutions,  the  form  of  the  algebraic  Reynolds  stress  and  dimensionality  tensors  needs 
to  be  determined  for  the  same  flows.  In  all  cases  the  direction  of  the  mean  and/or  frame 
rotation  rate  vectors  is  chosen  to  be  aligned  with  the  £3  direction.  Also,  in  all  cases  the 
eddy-axis  tensor  approaches  the  one  component  (1C)  point,  a2  =  1  in  the  limit  of  infinite 
total  strain,  when  the  exact  solutions  apply.  Then  in  principal  coordinates  of  the  eddy-axis 
tensor, 

an  =  1,  633  -  1.  (3.42) 
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Taking  =  — 1,  results  in  the  following  ASBM  predictions: 

dn=0,  di2=0,  d22  =  |(l-x)5  rf33  =  |(l  +  x)j  (3.43) 

rn  =  <t>,  n 2  =  -|7(1  +  x)>  ^22  =  J(l“0)(l  +  x)i  r33  =  §(1  -  0)(1  -  x)-  (3.44) 

The  model  for  7,  Equation  (3.26),  can  be  better  understood  now.  It  is  chosen  based  on 
the  observation  that  in  the  limit  of  infinite  total  strain,  in  RDT,  when  7*12  does  not  vanish, 
it  correlates  well  with  ^7*11^22  (see  Equations  3.33,  3.38,  and  3.41).  f3  is  then  introduced  as 
a  parameter  in 

ri  2  =  /Vrnr22.  (3.45) 

Using  the  ASBM  representation  (Equations  3.44)  in  both  sides  yields  the  form  chosen  for 
7, 

37(1  +  x)  =  Py/<t>h(1-4>)(1  +  x), 

7  =  /V20(W)/(1  +  X).  (3.26) 

Therefore,  instead  of  modeling  7  directly,  the  expression  above  is  used  in  the  hope  that 
modeling  (3  is  much  simpler  (in  fact,  in  the  a2  =  1  plane  (3  =  1  or  @  0,  depending  on 

whether  the  turbulent  kinetic  energy  is  expected  to  increase  or  not). 

3.9.4  Structure  Scalar  Models 

The  structure  parameters  are  defined  in  a  three-dimensional  space,  in  terms  of  r)mi  77/, 
and  a2,  which  were  introduced  in  Equations  3.27,  3.28,  and  3.29.  Comparing  the  2D  RDT 
exact  solutions  for  the  Reynolds  stresses  and  the  turbulence  structure  from  Section  3.9.2, 
with  the  corresponding  ASBM  solutions  from  Section  3.9.3,  gives  models  for  the  structure 
parameters  in  the  limit  of  infinite  distortion. 

The  structure  parameters  are  first  defined  in  a  generic  a2  plane,  along  the  rjf  and  r\m 
directions,  trying  to  match  the  2D  RDT  solutions  when  applicable.  The  structure  parame¬ 
ters  are  then  specifically  sensitized  to  the  degree  of  anisotropy  of  the  turbulence,  measured 
along  the  a2  direction.  In  this  final  step  an  adjustment  is  made  to  ensure  the  structure 
parameters  meet  minimal  function  continuity  requirements  (found  by  requiring  continuous 
gradients  of  the  turbulent  shear  stress  in  shear  flows,  e.g.  a  channel  flow). 

In  the  generic  a2  plane,  the  structure  parameters  are  defined  along  the  mean-shear  line, 
r)m  =  1,  and  along  the  plane-strain  line,  rjm  =  0,  for  a  generic  a 2  plane.  The  structure 
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parameters  are  then  interpolated  or  extrapolated  in  the  same  a2  plane,  depending  on  the 
flow  location  in  this  a 2  plane,  specified  in  terms  of  rjf  and  rjm. 

In  the  following,  the  subscripts  “0”  and  “1”  applied  to  0,  /?,  and  x,  refer  to  values  along 
the  lines  r]m  =  0  and  rfm  =  1,  respectively.  The  superscript  is  used  to  denote  values  on 
the  a2  plane  where  rjm  =  0  and  r\m  —  1  were  evaluated.  The  final  symbol  for  the  structure 
parameters  is  devoid  of  indices. 


Generic  a2  Plane,  rim  =  1  (Mean  Flow  Corresponds  to  Shear). 


Comparing  the  RDT  solution,  Table  3.1,  to  the  ASBM  solution,  Equations  (3.43)  and 
(3.44),  leads  to  the  forms  below.  Recall  that  the  subscript  “1”  indicates  values  computed 
along  the  line  rjm  —  1. 

If  rjf  <  0, 


,  (3.46) 


If  0  <  T)f  <  1, 

*,  =  (!-„),  ft-1,  +  (3-47) 


If  7]f  >  1, 


1 


,  73  .ox 

*'  “  1  (1  —  °2)  +  (*7/  —  1)’  <3'48) 


where 

6o  =  1.0  61-IOO  62-O.8  (3.49) 


were  optimized  for  channel  flow  with  spanwise  frame  rotation  (aligned  with  mean  flow 
rotation).  Notes: 


•  0i  (rjf  — >  ±00)  1/3  instead  of  1/2  (suggested  by  Table  3.1)  so  that  7*12  — ►  0  when 

rjf  — *  ±00,  regardless  of  a^-.  In  RDT,  an  — *  1  and  ri2  — >  0  regardless  of  0.  However, 
in  more  general  flows  it  is  possible  to  have  rjf  — ►  ±00  with  finite  deformations  (for 
example  near  the  zero  velocity  gradient  position  in  a  rotating  channel  flow).  In  this 
case,  coupling  0  — ►  1/3  with  /3,%  — >  0  implies  ri2  — >  0  regardless  of 
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•  factors  of  (1  -  a2)  enforce  0  — >  0  when  appropriate  (RDT  limit); 

•  factors  of  (rj/^/a2  -  1/3  oc  & r)  are  used  to  decrease  0\,  especially  at  small  rjf  <  0  or 
rjf  >  1,  for  larger  frame  rotation  rates  (consider  Tr  oc  yja2  -  1/3  and  rjf  oc  fi^r/IY). 
A  smaller  0  is  needed  to  get  the  large  asymmetries  observed  in  spanwise  rotating  chan¬ 
nel  flows.  Originally  a  large  constant  was  used  instead  of  this  form,  but  while  that 
gave  a  good  asymmetry  for  the  high  rotation  numbers,  it  also  displayed  a  large  asym¬ 
metry  at  small  rotation  numbers.  Hence  a  function  that  vanishes  at  small  rotations 
was  chosen. 

•  for  simplicity  xi  is  defined  in  terms  of  0\  (except  where  indicated  otherwise). 


Generic  a2  Plane,  r?m  =  0  (Mean  Flow  Corresponds  to  Plane  Strain). 


Comparing  the  RDT  polynomial  approximations,  Equations  (3.36)  and  (3.37),  to  the 
ASBM  solution,  Equations  (3.43)  and  (3.44),  leads  to  the  following  (recall  that  the  subscript 
“0”  is  used  to  indicate  values  computed  along  the  line  r)m  =  0): 

If  rtf  <  v^/4, 


(3.50) 

(3.51) 

(3.52) 


else,  if  T]f  >  \/3/4, 


0o  = 


(Vf-W 3/4), 
(1-a2) 


-l 


Xo  =  -0o, 

0o  =  (1  —  A>)/3, 


(3.53) 

(3.54) 

(3.55) 


where  63  =  1.0.  Equations  (3.53-3.55)  are  arbitrary  functional  forms.  The  objective  is  to 
retrieve  <j> 0  =  1/3  and  Xo  =  00  =  0  for  large  rjf  and/or  a2  — » 1.  The  2D-RDT  plane  strain 
solution  actually  applies  for  all  rjf,  predicting  0  —  0  and  x  =  — 1  for  ??/  >  \/3/4.  So  both 
the  functional  shapes  and  the  threshold  of  rjf  after  which  they  are  applied  are  arbitrary. 
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Generic  a2  Plane,  Interpolation  Between  rjm  =  0  and  =  1 


As  mentioned  earlier,  there  is  a  trapezoidal  region  in  the  (r)m,Vfia2  =  1)  plane  where 
the  turbulence  kinetic  energy  grows  for  two-dimensional  turbulence.  There  is  an  exact 
solution  for  the  Reynolds  stresses  and  for  the  structure  dimensionalities  (and  in  turn  for  the 
structure  parameters)  for  each  point  in  this  parallelogram,  which  however  is  not  in  a  simple 
algebraic  form.  It  involves  integrals  of  hypergeometric  functions.  To  circumvent  specifying 
values  at  each  point,  values  are  specified  along  the  edges  of  the  parallelogram,  and  then 
interpolated  inside.  The  interpolation  actually  is  done  for  every  point  between  the  lines 
r]m  =  0  and  T}m  =  1,  it  just  relates  to  an  exact  solution  inside  the  quadrangle;  outside  of 
it  the  exact  solution  predicts  decay  of  turbulence,  so  all  that  is  required  is  something  that 
drives  the  turbulence  to  zero,  notably  /?  — ►  0. 

The  interpolation  function  is  selected  based  on  the  simpler  case  of  combined  mean  strain 
and  mean  rotation  in  a  fixed  frame,  rjf  =  0,  for  which  an  exact  closed- form  RDT-solution 
exists  for  the  turbulence  structure.  The  variation  of  the  exact  solution  between  rjm  =  0  and 
rym  —  1,  for  other  values  of  77/  is  not  expected  to  differ  substantially.  Other  interpolation 
schemes  could  be  tried  in  the  future  if  deemed  necessary. 

Comparing  the  2D  RDT  solution  for  combined  mean  strain  and  mean  rotation  in  a  fixed 
frame,  Equations  (3.40a-f),  to  the  ASBM  solution,  Equations  (3.43)  and  (3.44),  gives  the 
variation  for  the  structure  parameters  along  0  <  r)m  <  1,  when  rjf  —  0, 


0  ~ 
X  = 


1(1  —  In  2)77^7(1  —  77^)  _  2.25 

1  +  1(1 -In 2^/(1-^)  ~  m 
1  +  f  (8  +  7T2  -  16 In 4  +  31n24)^/(l  -  rfa) 

1  +  I  (1  —  In  2)j?m/ (1  ~  Vm) 


0-185 r,2m25. 


(3.56) 

(3.57) 


Notice  the  discontinuity  between  the  value  of  x  for  shear  on  a  stationary  frame,  \  —  1/2, 
and  its  limiting  value  for  combined  strain  and  rotation,  \  ~  0.185.  On  this  work  the  former 
value  is  used  in  the  shear  limit.  The  approximate  profiles  are  used  nevertheless  as  guides 
to  the  interpolation  formulas. 

To  simplify  the  interpolation,  the  trapezoid  where  there  is  turbulence  growth  is  trans¬ 
formed  into  a  rectangle  with  the  transformation 


V*  =  -v m  +  [4/V3  +  (2  -  4/Vs)r]m}rif. 


(3.58) 


The  structure  parameters  are  then  interpolated  along  lines  of  constant  77*  (in  the  region 
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7]m  <  1)  according  to 


4>*  =  Mv*)  +  [Mv*)  - 

(3.59) 

F  =  Poiv*)  +  [Pl(V*)  ~  Po(V*))Vm, 

(3.60) 

X*  =  xoiv*)  +  ixiiv*)  -  Xo (v*)\vt- 

(3.61) 

Note  that  the  combined  mean  strain  and  rotation  solution  is  abandoned  in  favor  of  a  profile 
that  reaches  x  =  1/2  smoothly  at  the  shear  point. 


Generic  a2  Plane,  Extrapolation  for  7]m  >  1 


There  is  no  exact  steady  solution  for  the  region  rjm  >  1.  In  this  region  the  turbulence 
structure  probably  reaches  a  limiting  state  where  it  oscillates  periodically  (see  Blaisdell 
and  Shariff  1996).  On  the  other  hand  the  algebraic  model  responds  to  the  amount  of 
deformation  applied.  If  the  normalized  deformation  rate  is  steady,  the  response  of  the 
model  is  also  steady,  so  it  cannot  predict  a  limit  cycle.  Nevertheless  rotation  is  expected  to 
have  a  decorrelating  effect  on  the  turbulence,  so  the  structure  parameters  are  modeled  to 
reflect  this.  In  particular,  as  the  mean  rotation  grows  for  a  fixed  strain  — >  oo),  0  and 

X  are  driven  to  zero,  while  <£  is  driven  to  1/3  (for  consistency  with  the  case  of  infinite  r]f  in 
the  case  r)m  =  1  investigated  previously).  Hence  the  structure  parameters  are  extrapolated 
from  the  rjm  >  1  line  according  to 


<£*  =  1/3  + 


(&(*?/)  - 1/3) 

l  +  fom-l  )/(l-a2)’ 


Mm) 

l  +  fam-l  )/(l-a2)’ 
Xifa/) 

1  +  (Vm.  —  1)/(1  —  a2) 


These  are  also  arbitrary  functions. 


(3.62) 

(3.63) 

(3.64) 


a2  Direction  (Anisotropy  Characterization),  and  Function  Continuity 

So  far  the  structure  parameters  have  been  set  on  a  given  a2  plane.  For  <£i  this  is  the 
a2  =  1  plane,  except  for  a  possible  dependence  of  rjm  and  rjf  on  the  state  of  the  eddy-axis 
tensor.  The  other  scalars,  0\  and  Xh  have  been  sensitized  to  a2  in  some  form,  depending 
on  the  position  of  the  flow  in  the  (7?m,r?/)  plane.  At  this  point  there  are  two  final  conditions 
to  be  met: 
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(i)  in  isotropic  (undeformed)  turbulence  with  a 2  =  1/3,  the  structure  is  expected  to  con¬ 
sist  of  axisymmetric,  (x  —  0),  vortical,  {(j)  =  0),  eddies; 

(ii)  the  structure  parameters  </>,  7,  and  x>  need  to  be  Cl  continuous  functions  of  rjf , 
for  all  rjf.  When  rjf  — ►  rjm,  f3  needs  to  be  only  CO  continuous  with  respect  to  r?/, 
while  <j>  needs  to  be  0((rjf  —  r/m)3)  as  rjf  — ►  rjm.  These  requirements  arrive  because  of 
the  requirement  that  the  Reynolds  stresses  be  Cl  continuous,  such  that  their  spatial 
gradient  in  the  momentum  equation  be  CO  continuous. 


These  conditions  are  explored  in  more  detail  at  the  end  of  this  section,  looking  at  the 
particular  case  of  shear  in  a  rotating  frame. 

The  final  model  for  the  structure  parameters,  satisfying  these  conditions,  is  (for  near- wall 
regions  see  Equations  3.80  and  3.81) 

(Vm  -  Vf)2  \(  IVm-VflVW1  -  5) 


<P  =  4>* 

(3  =  13*, 


a2)2) 


(Vm  -  Vf)2  +  (1  -  a2)2  J  \jr7m  -  ?7/|i/f(a2  -  |)  +  Po(l  ~  a 2), 

J3/  2  nr 

X  =  ^  [2{a  -b)\  > 


(3.65) 

(3.66) 

(3.67) 


where 


p0  =  1.3  Pl  =  1.0  (3.68) 

were  optimized  along  with  ao  (Equation  3.15)  in  order  to  satisfy  the  key  condition  on  a 
boundary  layer  in  zero  pressure  gradient  observed  experimentally, 

=  —0.15,  simultaneously  with  —  =  —  4r12— ^  =  1,  (3.69) 

qz  e  e 

where  X\  and  #2  are  respectively  the  streamwise  and  wall-normal  directions.  Further  guid¬ 
ance  for  the  model  constants  was  provided  by  the  normal  stresses  (rn  =  0.53,  =  0.19, 

and  7-33  =  0.28  are  representative  values  in  the  log  region  of  the  DNS  of  Kim  (1992),  for  a 
channel  flow  at  Rer  —  385). 

Notes: 


•  the  factor  [3/2(a2  — 1/3)]  is  an  interpolant  between  the  isotropic  state,  when  it  is  zero, 
and  the  1C  condition  for  aij  (a2  =  1),  when  it  is  one; 

•  the  factors  of  (1  —  a2)  in  </>  are  such  that  the  RDT  form  is  recovered  for  a2  =  1,  while 
the  factors  of  (r)m  —rjf)  give  the  proper  degree  of  continuity  required  of  (j)  such  that 
item  (ii)  above  is  satisfied. 
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To  elaborate  in  the  functional  form  used  in  Equations  3.65  and  3.67  consider  the  case 
of  shear  in  a  spanwise-rotating  frame.  The  mean  deformation  is  given  by  Equation  (3.32), 
the  flattening  tensor  is  b,j  =  6^5 and  the  eddy-axis  tensor  becomes 

0  N 
0 

a  33 

The  normalized  Reynolds  stress  tensor  then  is  (recall  ^33  =  1  —  rn  —  ^2) 

n  1  =  (1  -  4>) |(1  -  an)  +  fan  +  (1  -  <t>)x  [3(1  -  °33)  -  £(!  +  a33)au] 

+  (-7^3/fiT)a2l[l  +  x(l  -  O33)], 

^22  =  (1  —  <^)|(1  —  «22)  +  4>a22  +  (1  —  <f>)X  [§(1  —  a33)  —  |(1  +  a33)a22] 

+  (7^3  /^T)a12[l  +  X(1  —  “33)]) 
ri2  =  “(I  -  0) sa12  +  <Ml2  +  (1  -  <t>)X  [“M1  +  «33)ai2] 

+  (7fiI7ftT)(an  -  022)1  [1  +  X(1  -  033)]- 

A  first  observation  is  that  in  a  fixed  frame  of  reference  with  positive  T  =  G12,  the  ASBM 
(Equations  3.15  and  3.16)  predicts  <122  =  <*33  <  ail  and  <212  >  0.  Additionally,  when  Tr  0, 
the  ASBM  predicts  an  -  <222  =  0((IY)2)  and  ai2  =  O(rr).  Inspection  of  the  constitutive 
relation,  Equation  (3.71c)  reveals  that  if  <j>  >  1/3  as  Tr  — ►  0,  then  7*12  >  0,  i.e.  a  negative 
eddy  viscosity  in  an  equivalent  Boussinesq  approximation,  which  should  not  happen  in  this 
flow  (homogeneous  shear  in  a  fixed  frame  of  reference).  Enforcing  0  — >  0  when  Tr  — >  0 
prevents  it  from  happening. 

The  continuity  requirements  of  the  structure  parameters  were  analyzed  in  detail  for 
the  case  of  shear  flow  in  a  rotating  frame,  corresponding  to  all  the  wall-bounded  flows 
investigated  in  the  present  work.  Particular  attention  is  given  to  the  shear  stress  7*12. 
For  a  given  frame  rotation  rate,  condition  (ii)  above  requires  continuity  of  dr  12 / dr) f  (rjf). 
Inspection  of  equation  (3.71c)  shows  that  this  in  turn  requires  (£(77/),  x(Vf)i  and  /) to 
Cl  continuous.  Cl  continuity  of  7  requires  in  turn  (3(r)f)  to  be  Cl  continuous  (if  0/0)  and 
y/4>  to  be  Cl  continuous.  Along  the  shear  line  (77m=l)  Equations  (3.46),  (3.47),  and  (3.48) 
were  made  piecewise  continuous.  The  boundary  at  r\f  =  0  is  not  an  issue,  because  no  two 
adjacent  points  in  a  channel  flow  cross  this  boundary.  In  fact,  77/  — ►  ±0  as  either  wall  of  a 
channel  is  approached.  The  boundary  at  rjf  =  1  definitely  is  crossed  somewhere  in  the  flow. 


(3.70) 

(3.71a) 

(3.71b) 

(3.71c) 


an  a  12 
ai2  022 
0  0 
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Because  <j>  —  0  at  this  point,  f3  need  only  be  CO  continuous.  It  can  be  shown  by  expanding 
0  in  a  Taylor  series  that  the  Cl  continuity  requirement  on  a/0  demands  0  =  0((rjf  —  l)3) 
as  rjf  — »  1.  Equation  (3.65)  was  constructed  to  achieve  just  that. 

3.10  Wall  Blockage 

As  a  no-slip  wall  is  approached,  the  tangential  components  of  the  velocity  are  driven 
to  zero  through  the  action  of  viscous  forces.  Additionally,  the  velocity  vector  is  reoriented 
into  planes  parallel  to  the  wall  through  an  inviscid  mechanism  (wall  blocking)  which  acts 
over  distances  far  larger  than  the  viscous  length  scale.  The  end  result  is  that  the  velocity 
component  normal  to  the  wall  is  driven  to  zero  faster  than  the  tangential  components. 

Here  wall  blocking  is  modeled  through  the  use  of  a  modified  projection  operator,  applied 
to  the  tensor  (and/or  vector)  to  be  blocked,  as  in  Reynolds  et  al  (2000).  The  strength  of 
the  projection  depends  on  the  proximity  to  the  nearest  wall,  as  measured  by  a  scalar. 

In  the  structure-based  model  it  is  postulated  that  at  the  wall  all  eddies  shall  he  in  the 
plane  of  the  wall.  Wall-blocking  is  then  introduced  to  reorient  the  eddies  into  planes  parallel 
to  the  wall.  The  structure  parameters  are  also  sensitized  to  wall  blocking,  such  that  the 
velocity  components  induced  by  the  eddies  are  consistent  with  the  asymptotic  behavior  or 
the  turbulence  velocities  near  walls. 

If  all  eddies  are  parallel  to  the  wall,  the  jetal  mode  of  motion  represents  motion  strictly 
parallel  to  the  wall.  The  vortical  mode  however  typically  involves  a  velocity  perpendicular 
to  the  wall.  To  eliminate  this  contribution,  two  different  approaches  can  be  proposed.  The 
first  regards  flattening  of  the  vortical  motion  around  the  eddies.  If  vortical  motion  were 
to  happen  on  a  fully  flattened  eddy  (resembling  a  sheet  on  the  plane  of  the  wall),  then  no 
normal  velocity  would  be  generated.  Alternatively,  one  can  postulate  the  absence  of  vortical 
motion  at  the  wall,  allowing  only  jetal  motion.  This  also  prevents  the  existence  of  a  normal 
velocity  at  the  wall. 

It  would  seem  that  a  less  intrusive  form  of  blocking  would  enforce  flattening  of  the 
vortical  part  of  the  motion  only.  However,  this  would  require  some  form  of  conditional 
sampling  and  flattening,  which  increases  the  complexity  of  modeling.  Simpler  alternatives 
are  (i)  fully  flattening  all  eddies  on  the  plane  of  the  wall  or  (ii)  enforcing  fully  jetal  eddies. 

Full  flattening  is  not  attractive,  because  of  its  secondary  effect  of  causing  cfe  —*  1  (taking 
X2  to  be  the  direction  normal  to  the  wall).  The  fully  jetal  alternative  is  more  palatable, 
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as  least  in  the  presence  of  shear.  Shear  flows  in  the  presence  of  walls  display  mostly  jetal 
motion  anyway,  so  enforcing  0=1  does  not  change  the  character  of  the  eddy  representation 
substantially.  However,  applying  the  model  in  shear-free  wall-bounded  flows  may  require  a 
reevaluation  of  this  approach. 

Blockage:  Modeling  Approach 

Following  Reynolds  et  al  (2000),  the  homogeneous  eddy-axis  tensor,  a-,  is  computed 
based  on  the  homogeneous  algebraic  procedure,  Equations  (3.15)  and  (3.16)  (note  that  the 
superscript  “h”  has  been  added  in  the  current  section).  It  is  then  partially  projected  onto 
planes  parallel  to  the  wall, 

aij  =  HikHjidki,  (3.72) 

where  the  superscript  h  indicates  an  homogeneous,  unblocked,  quantity.  Similarly  to  the 
regular  projection  operator  (see  Strang  1993),  the  partial-projection  operator  is  defined  as 

B-ik  =-—(Sik  —  Bik)i  =  1  —  (2  —  £/efc)a^nI?nm,  (3.73) 

L>a 

where  D2  is  such  that  the  trace  of  remains  unity.  The  blockage  tensor  Bij  gives  the 
strength  and  the  direction  of  the  projection.  If  the  wall-normal  direction  is  then  B22 
is  the  sole  nonzero  component,  and  varies  between  0  (no  blocking)  far  enough  from  the 
wall,  to  1  (full  blocking)  at  the  wall.  If  B22  =  1  —  0(y)  as  y  — ►  0  (approaching  the  wall), 
then  the  “2” -component  of  the  blocked  tensor  is  dampened  as  0(y).  Thus  ai2  =  0(y)  and 
a2 2  —  0(y2),  while  the  trace  an  =  1  is  preserved. 

Inspired  by  Durbin  and  coworkers  (Durbin  and  Pettersson  Reif  2002;  Parneix  et  al 
1998),  B^  is  computed  by  a  method  involving  an  elliptic  relaxation  equation, 


d2$ 

L—  =  23, 

(3.74) 

dxkdxk  ~ 

V 

at 

solid  boundaries, 

(3.75) 

.£-« 

oxn 

at 

open  boundaries, 

(3.76) 

where  xn  is  the  direction  normal  to  the  boundary.  B^  is  computed  by 

Bv  =  7Tir-*  if  ****>  0. 


(3.77) 


36 


If  all  gradients  of  $  vanish,  B^  is  computed  from  an  average  over  surrounding  points. 
Expanding  on  a  Taylor  series, 

Bij  =  if  9Jt9ik  =  0-  (3-78) 

The  length  scale  in  the  elliptic  relaxation  equation,  (3.74),  is  modeled  as  a  constant:  a 
first  order  model.  In  a  parallel  channel  flow,  this  implies  that  B22  becomes  a  hyperbolic 
cosine  function  of  the  distance  to  the  wall.  However  in  complex  flows  uT  (used  for  nor¬ 
malizing  L  in  Equation  3.74)  may  be  ill  defined.  For  instance,  near  the  separation  and 
reattachment  points  of  a  complex  flow  uT  vanishes.  Additionally,  in  a  complex  flow  region 
( e.g .  near  a  corner)  it  may  not  be  clear  what  should  be  used  for  uT.  Therefore,  in  general 
flows  a  more  elaborate  model  for  T,  in  terms  of  local  turbulent  (fc3/2/e)  and  Kolmogorov 
((i^/e)1/4)  length  scales,  should  be  adopted. 

The  specific  value  for  L  was  chosen  to  reproduce  the  intercept  C  in  the  logarithmic  law 
of  the  wall  for  the  mean  velocity, 

U/uT  =  (1  /«)  log {yur/v)  +  C 

where  U  is  the  local  mean  velocity,  u2  is  the  wall  shear  stress  divided  by  the  density,  and 
y  is  the  distance  to  the  wall.  The  blocking  modification  reduces  the  Reynolds  shear  stress 
near  the  wall,  with  the  amount  of  reduction  depending  on  the  value  of  L.  L  then  has  a 
strong  influence  over  the  intercept  C  in  the  log  law.  Away  from  the  wall,  the  wall  blocking 
modification  dies  out,  so  it  does  not  affect  the  slope  of  the  log  law. 

In  Reynolds  et  al  (2000)  the  projection  operation  was  also  performed  on  the  homo¬ 
geneous  Reynolds  stresses.  Here  this  approach  is  abandoned,  and  instead  the  structure 
and  the  Reynolds  stresses  are  computed  based  on  blocked  ay,  0,  and  7^.  This  generates 
realizable  tensors  even  in  presence  of  wall  blocking. 

Turning  back  to  the  channel  flow,  with  B22  =  1  —  0(y),  a\2  =  0(y ),  and  022  =  0(y 2), 
the  blocked  normalized  Reynolds-stress  components  (Equations  3.71b, c)  become 

r22  =  (1  -  <p)i  +  (1  -  4>)x  [5(1  -  033)]  +  ml/SF)an[l  +  x(l  -  033)]  +  0(y2),  (3.79a) 

m  =  (7^3  /nTVn  2 11  +  X(1  -  033)]  +0(y).  (3.79b) 

The  proper  asymptotic  behavior  of  the  Reynolds  stresses,  rn  oc  0(y)  and  7*22  oc  0(y2),  is 
recovered  if  the  jetal  parameter  behaves  as  0  =  1  — 0(y2),  and  the  helix  parameter  7  behaves 
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as  7  =  0(y).  To  achieve  that,  they  are  blocked  using 

<£  =  1  +  (<^  -  1)  (1  -  Bkkf ,  (3.80) 

7  =  7h(l-£u),  (3-81) 

where  <j)h  and  are  respectively  the  homogeneous  jetal  parameter  and  the  homogeneous 
helix  parameter.  The  trace,  is  used  instead  of  the  component  in  the  wall-normal 
direction,  Bnn  (no  summation)  where  xn  is  the  wall-normal  direction.  They  are  both 
1  —  0(y )  near  the  wall,  but  Bkk  =  $  is  simpler  because  it  does  not  require  the  specification 
of  a  wall-normal  direction  everywhere  in  the  flow.  Away  from  the  wall,  when  Bkk  —>  0, 
both  parameters  approach  their  homogeneous  values. 

A  consequence  of  this  approach  is  that  as  the  structure  tensors  are  derived  from  the 
eddy-axis  tensor,  realizability  is  automatically  satisfied  for  ry,  dy  and  /y.  Furthermore 
the  fundamental  constitutive  equation  relating  the  three  tensors  remains  satisfied,  as  the 
representation  of  the  vector  triad  (tf,  s,  n)  in  eddy-axis  coordinates  remains  unaltered. 

3.11  Model  Summary 

Given  a  mean  strain-rate  tensor,  mean  rotation-rate  tensor  and  frame-rot  at  ion-rate  ten¬ 
sor,  Siji  Clij,  and  fl{j  respectively,  with  the  total  rotation-rate  tensor  fiy  =  fiy  +  fl{j,  and 
a  time  scale  of  the  turbulence,  r,  compute: 

3.11.1  Homogeneous  Structure 
Strained  afj. 


■*IJ 


0%  + 


(Sjkakj  +  ^jkaki  ~  %^mnanm^ii)T 


3  00  +  2^1-2 
where  S*j  —  Sy  —  Skk^ij/ 3  is  the  traceless  strain-rate  tensor. 


Rotation  Operation,  Producing  the  Homogeneous  Oy. 


(3.15) 


aij  —  HikHjiah, 


(3.16) 
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where 


Hij  =  5ij  +  hi  — r==  +  h 2 


O2  ’ 

*Lpp 


where  Opp  =  Qpqttpq,  and 


hi  =  J2h2  -  tq/2, 


2~2Ji(l  + vT^r)  if  r  <  1 


2-2^j(l- v/l  =  l7r)  if  r  >  1 , 

Clpq£lqrSTp 

Sir>Snmamk 


Flattening  Tensor  Model 


b  (ni  +  cbn{)(uj  +  cbnfJ) 
ij  (fifc  +  cbnfk)(nk  +  cbn{)  ’ 

where  and  Q{  are  respectively  the  mean  vorticity  and  the  frame  vorticity  vectors. 


Helix  Vector  Model 


Ik  =  7 


1 20(1  —  </>) 


(3.25-3.26) 


Structure  Scalars:  Parameterization 


Vf  —  Vm  —  sign(X)y 


(Z  —  Q'pqQ'pqi 


(3.27-3.29) 


—  O'ijQik^kj >  ^  —  ^ij^ik^kj^  X  _  aijQ^Skj-  (3.30) 
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Structure  Scalars;  a2  Plane,  7jm  =  1 

If  r]f  <  0, 

01  =  111,  ft  =  - ^ - 

1  -  ‘"Tr^T  (‘ + 


’»  X  2^lj 


else,  if  0  <  77/  <  1, 

0i  =  (!-»?/),  A  =  1>  Xi  =  2  +  2  (*  “  1  +  M/Al  -  a2) )  ’ 
else,  if  7]f  >  1, 

,  ns- 1  *  1  ..  ,  (1  -  Pi)(vs  - 1) 

*"-1’  A  "(1-o!,+fe-1 

Structure  Scalars;  a2  Plane,  77™  =  0 

If  Vf  <  V3/4, 


Vl  - 1  v  /o  ^ 

(i-.»)+(w-i)-  <3-48) 


'(2i?/)2 

a%)2V1 

.  3/4 

\  3/4  /  J 

00  =  0.145 
0o  =  h 


Xo  =  -  0.342^|^ +  (1-0.342) 


(a  vsr\6' 

V  3/4  ;  j  ’ 


else,  if  77/  >  a/3/4, 


0o=  1  +  6»'(!?{i^W(°2-»  > 

Xo  =  -/?o, 

0o  =  (l-/3o)/3. 

Structure  Scalars;  a2  Plane,  Interpolation  Along  the  (97 miWf)  Directions 

If  rr]m  5: 

V*  =  ~~Vm  +  [4/ a/3  +  (2  —  4/ \/3 
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<t>*  =  Mv*)  +  IMv*)  -  Mv*)]vh, 
P*  =  Po(v*)  +  [Pi(v*)  -  Po(v*)}Vm> 
X*  =  Xo(v*)  +  lxi(v*)  -  Xo(r)*)]Vmi 


else,  if  rjm  >  1 


^-1/9  ,  (Mvf)-m 

*  -1/3+l  +  (r,m-l)/(l -a*)' 


P*  = 


X  = 


Pi(Vf) 


1  +  (»? m-  1)/(1  -  a2)  ’ 
Xi (Vf) 


l  +  (T]m-l)/(l-a2)' 

Structure  Scalars,  Final  Homogeneous  Expression 


(3.59) 

(3.60) 

(3.61) 


(3.62) 

(3.63) 

(3.64) 


* 

II 

(  (Vm  ~  Vf)2  \ 

(  \Vm-VfW\{a2  -\)  | 

(3.65) 

\(Vm  -  Vf)2  +  (1  -  a2)2) 

\\Vm  —  VfWKa2  ~  %)  +Po(l  ~  fl2) /  ’ 

II 

* 

(3.66) 

X  =  X* 

w-w- 

(3.67) 

3.11.2  Wall  Blockage 

In  what  follows  a  subscript  h  indicates  a  homogeneous  quantity,  as  computed  above. 
The  homogeneous  quantities  are  modified  for  near  wall  effects  as 

aij  =  HikHjia.%1,  (3.72) 

Hik  =-^-(. 5ik  -  Bik),  D2a  =  1  -  (2  -  Bkk)ahmnBnm,  (3.73) 

Lfa 


<fi=l  +  (<j>h-l)(l-Bkk)2 
7  =  1h  (1  -  Bkk) , 


(3.80) 

(3.81) 
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where 


and 


Bij  —  \ 


if 

fjfe*  « 

=  %  i— =  23, 
oxkdxk  y 

4>  =  1  at  solid  boundaries, 

d$ 

&  n  =  — —  =0  at  open  boundaries, 

OXfi 

where  xn  is  the  direction  normal  to  the  boundary. 

Finally  the  structure  parameters  are  modified  as 

0=1  +  (<t>h  —  1)  (1  —  Bkk)2 , 

7  =  7h  (1  -  Bkk) , 


Note  that  x  is  not  modified  near  the  wall. 


(3.77-3.78) 


(3.74) 

(3.75) 

(3.76) 


(3.80) 

(3.81) 


Normalized  Structure  Tensors:  Reynolds  Stress,  Dimensionality,  and  Circulicity 


d ij  —  §[1  ~  X(1  —  ®nmbrnn)]fiij  (3.12) 

i  [-^  “h  Onm^mn)]aij  d"  X^ij  ~  X{ain^nj  4“  ajn^m)? 

rij  =  (1  -  <j>)l(8ij  —  dij)  +  <j)CLij  (3.13) 

+  (1  —  0)X  [§(1  —  Q>nmbmn)fiij  ““  £(1  +  Q>nmbmn)Q'ij  ~  b{j  +  &inbnj  4"  ajnbni\ 

+  (“7f2 J / fiT) (eiprapj  +  ejpr®pi){|[l  —  X(1  “  ®nm^mn)]4r  H”  X^kr  ~  X^knbnr} > 
fij  —  $ij  ~~  rij  ~~  dij.  (3.14) 

Model  Constants 

The  set  of  model  constants  used  with  the  model  are 

clq  =  1.6,  Cb  =  — 1.0,  bo  =  1.0,  b\  =  100,  62  =  0*8,  63  =  1.0,  po  —  1.3,  p\  =  1.0. 

(3.82) 
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3.12  Response  of  the  Homogeneous  ASBM  to  Various  Nor¬ 
malized  Deformation  Rates 

The  behavior  of  the  algebraic  model  equations  that  determine  the  turbulence  structure 
and  Reynolds  stresses  is  investigated  for  various  normalized  deformation  rates,  correspond¬ 
ing  to  different  mean  flows  and  different  ratios  of  the  turbulent  time  scale  to  the  mean 
flow  time  scale.  In  particular  the  model  response  to  frame  rotation  is  sought.  It  cannot 
be  overemphasized  that  the  results  from  this  section  do  not  correspond  to  a  time  inte¬ 
gration  of  the  turbulence  equations.  The  algebraic  equations  for  the  turbulence  structure 
(section  3.11.1)  return  a  turbulent  state  to  a  given  mean  normalized  deformation  rate, 
irrespective  of  any  time  history.  Results  here  display  the  response  of  the  algebraic  equa¬ 
tions  (the  algebraic  model  for  the  structure)  to  different  combinations  of  mean  normalized 
deformation  rates. 

3.12.1  Irrotational  Flows 

Flows  with  no  mean  and  frame  rotation  ( e.g .  plane  strain  or  axisymmetric  contrac¬ 
tion/expansion)  correspond  to  the  (0,0)  point  in  the  plane,  and  the  algebraic 

expressions  for  the  structure  scalars  return  (j>  =  x  =  0  (implying  7*.  =  0)  irrespective  of  the 
normalized  strain  rate  (which  is  measured  here  by  a2).  The  constitutive  relations  (3.12-3.14) 
then  return 


dij  —  rij  —  \  (Sij  —  a,ij ),  fij  —  Oij,  (3.83) 

as  predicted  in  Kassinos  and  Reynolds  (1994). 

3.12.2  Axisymmetric  Flows  in  a  Rotating  Frame 


The  mean  deformation  is  given  by 


(  r 

0 

0  ^ 

f° 

0 

0  ^ 

( 0 

0 

0  \ 

c?.  — 

DIJ  — 

0 

-r/2 

0 

,  Qij  — 

0 

0 

0 

5  &{j  — 

0 
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& 

1° 

0 

-r/2 ) 

1° 

0 

1° 

-Qf 

0  J 

with  T  >  0  for  axisymmetric  contraction,  and  T  <  0  for  axisymmetric  expansion. 
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Figure  3.1  displays  b ij  =  —  Sij/ 3,  the  normalized  anisotropy  of  the  Reynolds  stress 
tensor,  versus  normalized  strain  rate,  with  and  without  imposed  frame  rotation  (Qf  =  2|F|). 
Figure  3.1(a)  explores  fyj  for  the  axisymmetric  contraction  case.  In  the  fixed  frame  the 
turbulence  is  completely  vortical,  and  as  the  eddy-axis  tensor  aligns  with  the  direction  of 
positive  strain,  the  axial  component  of  the  Reynolds  stress  decays  rapidly,  Figure  3.1(a). 
In  presence  of  frame  rotation  the  jetal  mode  is  activated,  and  the  axial  component  of  the 
Reynolds  normal  stress  does  not  decay  as  fast  initially.  Eventually,  for  very  large  Fr,  the 
eddy-axis  tensor  becomes  aligned  with  the  direction  of  rotation  and  the  turbulence  becomes 
materially  indifferent  to  the  rotation,  and  the  Reynolds  stress  components  become  the  same 
as  without  rotation. 

Figure  3.1(b)  shows  bij  for  the  axisymmetric  expansion  case.  Again  the  eddy-axis  tensor 
aligns  with  the  direction  of  positive  strain,  now  orthogonal  to  the  axis  of  symmetry.  In  the 
fixed  frame  case  (no  jetal  mode),  in  the  limit  of  very  large  normalized  deformation  rates, 
the  Reynolds  stresses  become  one-component.  In  presence  of  frame  rotation,  the  jetal  and 
helical  modes  are  activated,  and  the  anisotropy  of  the  stresses  does  not  grow  so  strong.  In 
this  case  the  eddy-axis  tensor  never  becomes  aligned  with  the  direction  of  the  frame  rotation 
(a22  =  U33  — ►  1/2  as  Tr  — >  00),  so  the  turbulence  is  affected  by  the  rotation  at  all  levels  of 
Tr. 


3.12.3  Shear  Flow  in  a  Spanwise  Rotating  Frame 

The  mean  deformation  is  given  by 


S*.= 


/  0  r  0  ^ 


r  0  0 

000 


— 


(  0  r  0  ^ 


-r  0  0 

000 


Qf.= 


/  0  0  \ 

&  0  0 

000 


(3.85) 


Figure  3.2  shows  the  normalized  Reynolds  stress  versus  mean  normalized  strain  rate, 
for  several  ratios  of  frame  rotation  to  mean  rotation.  The  amount  of  frame  rotation  in  the 
different  figures  is  quantified  by  rjf  =  — n{2/Q{2,  the  ratio  of  the  (counter) -rot  at  ion  of  the 
frame  to  the  rotation  rate  of  the  mean  flow  (shear).  In  all  cases  the  eddy-axis  tensor  reacts  in 
precisely  the  same  way  to  the  different  normalized  strain  rates  because  the  eddy-axis  tensor 
is  computed  in  the  rotating  frame;  as  it  is  affected  only  kinematically  by  frame  rotation,  it 
is  insensitive  to  the  amount  of  frame  rotation.  However  the  stress  states  are  very  different 
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for  the  different  rjf  because  the  stresses  are  dynamically  affected  by  the  frame  rotation.  In 
the  ASBM  this  effect  comes  through  the  structure  parameters.  The  structure  parameters 
were  built  with  guidance  from  an  asymptotic  solution  to  RDT  of  the  same  problem  (shear 
with  spanwise  frame  rotation).  Hence  the  model  shows  the  same  kind  of  behavior  brought 
in  by  this  fitting. 

Note  that  for  zero  frame  rotation  the  r\\  component  is  largest,  but  as  rjf  increases  the 
other  normal  components  become  more  important.  Note  the  dramatic  effect  on  the  shear 
stress  ri2.  As  the  frame  counter-rotation  grows,  r\ 2  initially  grows  larger,  corresponding  to 
increased  turbulence  activity,  but  eventually  becomes  weaker  as  the  frame  counter-rotation 
overcomes  the  mean  rotation.  The  ratio  of  turbulence  production  to  turbulence  dissipation, 
P/e  =  ~4ri2rr  where  r  =  k/e,  is  directly  proportional  to  7*12,  and  once  the  ratio  becomes 
smaller  than  unity  turbulence  decay  is  observed. 

In  the  previous  model,  (Reynolds  et  al  2000),  r\2  would  change  sign  for  the  larger  values 
of  rjf  shown.  In  the  current  version  a  conscious  effort  was  made  to  postpone  this  to  larger 
rjf  and  normalized  strain  rates.  When  r\2  becomes  positive,  the  turbulence  production  rate 
becomes  negative!  More  discomforting  though  is  the  fact  that  this  amounts  to  a  negative 
effective  eddy  viscosity  (in  standard  two-equation  models  vt  oc  —  tWO^t)),  which  can  have 
a  destabilizing  effect  on  numerical  schemes.  In  particular  the  previous  ASBM  would  predict 
positive  7*12  for  any  normalized  strain  rate  if  rjf  was  large  enough.  Near  the  zero  velocity 
gradient  region,  in  spanwise  rotating  channel  flows,  rjf  — >  00,  and  the  model  solution  would 
predict  odd  oscillations  of  the  shear  stress  near  this  region  (as  rjf  decreases  r\2  would  change 
sign  again,  clearly  at  odds  with  the  DNS).  Great  efforts  were  taken  in  the  optimization  of 
the  structure  coefficients,  especially  /?,  so  that  this  situation  could  be  avoided. 
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For  the  particular  case  of  shear,  the  algebraic  eddy-axis  model  predicts  (Z22  =  <233.  Then 
using  the  constitutive  Equations  (3.13)  and  (3.12),  selected  Reynolds-stress  components  can 
be  explicitly  written  as 


n 3  =  -(1  -  4>)x  «22  &13  +  (7^3 /fiT)  «12  X  hz  (1  -  033)  (3.88) 

+  (7 012  [i(l  -  X)  +  5  X  (2&11  -  an  bn  +  a33  633)], 

^23  =  (1  -  <t>)x  °12  &13)  (3.89) 

^13  =  X  O22  ^13)  ^23  =  —X  O12  &13-  (3.90) 


The  sign  of  613  depends  on  the  sign  of  Cj.  Equation  (3.24)  becomes 


613  = 


qn{n3 

(0,  n{)2  +  (fi3)2' 


(3.91) 


The  case  of  shear  with  streamwise  rotation  is  very  challenging  and  exposes  some  diffi¬ 
culties  with  algebraic  models  in  general,  and  the  ASBM  in  particular.  The  Reynolds-stress 
tensor  becomes  fully  populated  as  the  presence  of  streamwise  rotation  activates  the  sec¬ 
ondary  shear  stresses  7*13  and  7*23.  It  is  not  clear  what  should  be  the  sign  of  these  secondary 
stresses.  Equations  (3.89)  and  (3.90)  show  that  the  signs  of  r23,  ^13,  and  d,2z  are  directly 
dependent  on  613  which  in  turn  depends  on  Q>.  When  the  frame  rotation  is  aligned  with 
the  mean  rotation,  the  sign  and  magnitude  of  Cb  are  immaterial. 

Shear  with  streamwise  rotation  does  not  benefit  from  the  guidance  of  the  asymptotic 
behavior  of  exact  solutions  (two-dimensional  turbulence  aligned  with  the  Xi  direction  as 
used  in  section  3.9.2  is  indifferent  to  rotation  around  the  aq-axis).  Nevertheless  some  guid¬ 
ance  can  be  obtained  by  comparing  the  relative  influence  of  streamwise  frame  rotation  on 
the  flow  by  means  of  the  Particle  Representation  Model  (PRM)  of  Kassinos  and  Reynolds 
(1994).  Specifically,  the  PRM  was  used  to  help  determine  the  proper  sign  of  the  secondary 
stresses  given  a  particular  frame  rotation  direction.  The  PRM  allows  the  time  integration 
of  turbulence  statistics  in  the  Rapid  Distortion  Limit  (RDT)  of  homogeneous  turbulence. 
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Figure  3.3(a)  compares  PRM  solutions  for  the  Reynolds  stress  tensor  in  the  no  frame  rota¬ 
tion  case  to  solutions  for  n23/^23  =  “0.5.  Figure  3.3(b)  gives  the  same  comparison  for  the 
dimensionality  tensor.  It  is  observed  that: 

•  the  effect  on  the  normal  stresses  and  the  shear  stress  ryx  is  noticeable  but  small.  Same 
for  the  dimensionality  tensor  components; 

•  ^13  >  0,  d2 3  <  0,  r23  >  0,  and  m  >  0  for  small  deformations. 

Figures  3.4(a)  and  (3.4(b))  compare  the  response  of  the  ASBM  for  the  Reynolds  stresses 
and  dimensionality  tensor  respectively,  for  different  normalized  rates  of  shear  in  a  fixed  frame 
and  in  streamwise-rotating  frame  with  fi23/^23  =  -0.5. 

A  direct  comparison  between  the  PRM  and  the  ASBM  results  is  not  possible,  as  they 
are  very  different  things.  The  PRM  represents  integration  in  time,  while  the  ASBM  stresses 
represent  solely  the  output  of  the  algebraic  equations  to  a  given  normalized  shear  rate  Tr. 
Nevertheless  relative  behavior  between  the  rotating  frame  case  to  the  fixed  frame  case  can 
be  investigated.  The  ASBM  results  were  obtained  with  C \>  —  —1.0.  This  gives  flattening 
&13  <  0  and  hence  secondary  Reynolds  stresses  and  dimensionality  components  with  sign 
consistent  with  the  PRM.  Larger  magnitudes  of  C&  do  not  affect  the  solution  substantially. 
Smaller  magnitudes  decrease  the  primary  shear  stress  further.  As  compared  to  the  PRM 
solution, 

•  the  principal  Reynolds  shear  stress,  ri2,  shows  excessive  sensitivity.  This  carries  over 
to  the  solution  of  streamwise  rotating  channel  flows,  which  will  display  too  much 
turbulence  suppression; 

•  the  secondary  stresses  are  much  smaller  then  what  is  attained  with  the  RDT  solution 
(PRM).  This  will  also  be  observed  later  in  the  streamwise  rotating  channel  flow. 

A  surprising  result  seen  in  the  dimensionality  curves  (PRM)  is  that  d22  >  c/33  for  the 
larger  times.  The  exact  solution  for  2D  RDT  predicts  —  3/4  and  d22  =  1/4.  This 
discrepancy  is  not  yet  understood. 

Part  of  the  deficiency  of  the  ASBM  should  be  attributed  to  a  carryover  effect  from  the 
spanwise  rotation  case.  The  structure  parameters  were  optimized  for  matching  the  observed 
suppression  of  turbulence  in  presence  of  spanwise  frame  rotation.  The  ASBM  however 
responds  to  rotation  through  scalars,  which  as  such  cannot  distinguish  the  direction  of  the 
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frame  rotation.  Thus  any  rotation  is  seen  as  stabilizing  rotation,  resulting  in  the  suppression 
of  the  shear  stress  in  the  ASBM.  Nevertheless  the  ASBM  can  be  less  susceptible  to  the 
streamwise  frame  rotation:  as  the  eddy-axis  tensor  aligns  with  the  streamwise  direction 
(due  to  shear),  it  becomes  less  sensitive  to  rotation  about  this  particular  direction  (material 
frame  indifference).  Standard  models  do  not  have  material  frame  indifference  built-in,  so 
they  become  more  susceptible  to  errors  because  of  the  optimization  for  the  spanwise  rotation 
case. 


3.13  Summary 

In  this  chapter  a  turbulent  field  is  represented  through  an  ensemble  of  hypothetical 
eddies,  with  assigned  properties  that  provide  a  measure  of  the  turbulence.  The  Reynolds 
stresses  and  the  structure  tensors  are  decomposed  onto  this  eddy  field,  resulting  in  consti¬ 
tutive  relations  for  the  Reynolds  stress  and  structure  tensors  in  terms  of  statistics  of  the 
eddy  field  —  the  eddy-axis  tensor,  the  eddy-flattening  tensor,  the  helix  vector,  and  scalars 
representative  of  the  eddy  structure.  The  definition  of  these  scalars  is  the  heart  of  the  chap¬ 
ter.  They  are  chosen  to  match  the  isotropic  turbulence  state,  and  selected  RDT  limiting 
states.  Intermediate  states  are  computed  through  blending  functions  chosen  with  care  to 
match  a  canonical  shear  state,  among  other  criteria. 

The  model  can  be  outlined  as: 

rij  —  r  ij  (flij  >  7,  x)?  dij  —  dij  (dij ,  bij ,  <j>, 7,  x)  ? 

where  (homogeneous  part  —  unsensitized  to  possible  presence  of  walls) 

-  a%(SijT,  Slijr),  )t ),  4>h,  7h,  *  =  /(Syr,  fiyr,  q{j,  afy). 

The  influence  of  a  wall  is  imparted  to  the  eddy  statistics  through 

dij  —  fa (dij, Bij),  0  —  g(<jfl,Bkk),  7  =  h( 7  ,•£?/.*.), 

where  is  a  local  measure  of  the  influence  of  a  wall. 

The  homogeneous  model  has  been  explored  for  diverse  normalized  mean-deformation- 
rates  (axisymmetric  deformations,  and  shear  in  presence  of  spanwise  and  streamwise  frame- 
rotation).  The  model  provided  realizable  and  physically  acceptable  states  for  the  deforma¬ 
tions  tested.  However,  in  the  streamwise  frame  rotation  mode,  the  model  displayed  excessive 
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sensitivity  to  the  frame  rotation.  This  is  most  likely  a  consequence  of  the  emphasis  in  using 
the  spanwise  frame  rotation  case  in  defining  the  eddy  structure  parameters,  to  the  detriment 
of  the  streamwise  case.  In  turn  this  is  blamed  on  the  fact  that  the  structure  parameters  are 
sensitive  to  rotation  through  scalars,  which  obviously  cannot  impart  a  direction  of  rotation 
in  its  final  representation. 
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(f)  t)f  =  -2 


Figure  3.2:  ASBM  stress  states  for  homogeneous  shear  in  a  frame  rotating  about  the  spanwise  axis: 
(a)  non-rotating  frame;  (b)-(e)  counter-rotating  frame;  (f)  co-rotating  frame  . 
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Figure  3.3:  Comparison  of  time  integration  of  turbulence  statistics  in  KDT  (PRM  of  Kassinos  and 
Reynolds  1994)  of  homogeneous  shear  in  a  fixed  frame  (solid  lines)  and  in  a  streamwise-rotating 
frame  with  =  -Qi2/2.  Dotted  lines,  dashed  line:  “13”  component.  Dot-dashed  line:  “13” 
component,  (a)  Reynolds  stresses  r^;  (b)  dimensionalities  dtj. 


rx  rt 


(a) 


(b) 


Figure  3.4:  Comparison  of  ASBM  response  to  different  normalized  shear-rates  in  a  fixed  frame 
(solid  lines)  and  in  a  streamwise-rotating  frame  with  q£3  =  —$li2/2.  Dotted  lines,  dashed  line:  “13” 
component.  Dot-dashed  line:  “13”  component:  (a)  Reynolds  stresses  (b)  dimensionalities  dij. 


Chapter  4 


Scale  Equations 


In  this  chapter  two  scale  equations  are  proposed  to  complement  the  algebraic  structure- 
based  turbulence  model  introduced  in  Chapter  3.  The  scale  equations  proposed  here  are 
the  extension,  for  inhomogeneous  wall-bounded  flows,  of  the  model  equations  suggested  by 
Reynolds,  Langer,  and  Kassinos  (2002)  (RLK02)  for  homogeneous  flows. 

The  extension  of  the  model  equations  is  such  that: 

•  it  returns  the  homogeneous  forms  of  RLK02  when  applicable; 

•  the  source  terms  in  the  evolution  equations  for  the  scales  have  proper  asymptotic 
behavior  when  approaching  a  solid  wall; 

•  the  model  scale  equations  are  consistent  with  the  existence  of  a  logarithmic  region  in 
a  wall-bounded  flow; 

•  the  model  scale  equations  should  not  show  any  spurious  behavior  near  the  edge  of  a 
shear  layer  (near  a  turbulent /non-turbulent  interface). 

RLK02  follow  Obukhov  (1941)  and  Singh  and  Lumley  (1971)  and  hypothesize  an  aver¬ 
aging  process,  denoted  by  brackets,  that  can  extract  the  large-scale  coherent  motions.  For 
example  this  could  be  a  POD  modal  decomposition  that  captures  the  energy-containing 
turbulence  and  the  associated  large-scale  vorticity.  Combined  with  the  time  (or  ensemble) 
average  (denoted  by  overbar),  gives  the  triple  decomposition  (see  Reynolds  and  Hussain 


1972)  for  the  velocity  field, 

Ui  ~Ui  +  Ui  +  u",  (4.1) 

(lli)  =  U{  +  Ui,  fai)  ^  ^i,  i^i)  =  ^5  (^*2) 

and  for  the  vorticity  field, 

Ui  =  Ui  +  +  uf( ,  (4.3) 

(uji)  =  Ui  +  uu  (u^  =  Ui,  (ui)  =0.  (4.4) 
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A  single  prime  was  used  in  Chapter  2  to  indicate  fluctuating  quantities  in  the  Reynolds 
decomposition.  Double  primes  are  used  here  to  distinguish  the  current  small-scale  fluctu¬ 
ations  from  the  former.  Then  an  overbar  over  a  variable  denotes  an  average  value,  a  tilde 
over  a  variable  denotes  a  large-scale  coherent  fluctuation,  and  a  double  prime  denotes  the 
remaining  incoherent  fluctuation.  It  is  assumed  that  the  averaging  operations  commute 
with  differentiation,  and  that  the  different  modes  are  not  correlated,  that  is, 

Uiu"  =  0,  UiUj  =  uiuj  =  0,  {uiUj)  =  Ui{ufj )  =  0.  (4.5) 

This  is  a  consequence  of  the  remaining  incoherent  scales  representing  uncorrelated  motion. 
Thus,  even  if  the  large-scale  motions  and  the  incoherent  motions  have  same  frequency  or 
wavelength,  they  do  not  have  any  phase  correlation,  so  that  the  average  product  vanishes. 


4.1  Overview 

In  preparation  for  choosing  appropriate  near-wall  scaling  factors  for  different  terms  of 
the  evolution  equations  for  the  turbulence  scalars,  a  summary  of  the  expected  asymptotic 
behavior  for  the  flow  variables  is  given  in  Section  4.2.  Sections  4.3  and  4.4  introduce  time 
and  length  scales,  respectively  used  in  the  normalization  of  the  mean  deformation  in  Chapter 
3,  and  in  scaling  source/sink  terms  in  the  scale  equations. 

Sections  4.5  and  4.6  introduce  the  evolution  equations  for  the  turbulent  kinetic  energy 
and  for  the  large  scale  turbulent  enstrophy  respectively,  and  looks  at  specific  model  choices 
for  the  different  terms  in  the  equations.  Section  4.7  briefly  points  proper  boundary  condi¬ 
tions  for  the  coupled  system  of  scale  equations. 

Section  4.8  then  looks  at  the  log  region  of  a  boundary  layer,  and  at  the  interface  between 
turbulent  and  non-turbulent  flows.  This  determines  conditions  to  be  satisfied  by  the  model 
constants,  such  that  the  model  conforms  to  accepted  behavior  in  these  regions. 

There  are  no  results  explored  in  this  chapter,  as  it  depends  on  the  specification  of  a  given 
turbulent  state,  which  is  done  in  this  work  using  the  algebraic  equations  from  Chapter  3. 
Full  model  results  (scale  equations  coupled  with  algebraic  structure  equations)  follow  in 
Chapter  5. 
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4.2  Near- Wall  Asymptotic  Behavior 

Taking  Xi  =  {x,y,z},  Ui  =  {u,v,w},  and  u,  =  {U(y),V  =  0,tV  =  0},  the  near-wall 
expansions  for  turbulent  fluctuations  and  statistics,  with  the  ^/-direction  as  the  wall-normal 
direction,  are 


u  =  hy+ciy2+  0(y3), 

(4.6a) 

v  =  c2y2+  0(y3), 

(4.6b) 

w  =  hy+C3y2+  0(y3), 

(4.6c) 

P  =  Po+  PiV+P2V2+  0(y3), 

(4.6d) 

U  =  Biy+Ciy2+0(y3), 

(4.6e) 

where  all  the  different  coefficients  are  functions  of  x  and  2.  Though  the  expressions  have 
been  written  only  for  the  large-scale  fluctuations,  they  apply  equally  to  the  small-scale 
fluctuations,  with  different  coefficients. 

It  follows  that  the  near- wall  expansion  for  the  fluctuating  vorticities  are 


Ux  =  wtV  -  vyZ  =  bs  +  2czy  +  0(y2), 

(4.7a) 

uv  -  uz  -  w,x  =  (bhz  -  bstX)y  +  0(y2), 

(4.7b) 

uz  =  V'X-utV=  -bi-2ay  +  0(y2). 

(4.7c) 

4.3  Time  Scale 

The  algebraic  equations  developed  in  Chapter  3  respond  to  the  mean  deformation  nor¬ 
malized  with  a  turbulence  time  scale,  r  (strictly  r  is  needed  only  for  normalizing  the  mean 
strain  rate  in  Equation  3.15  —  it  cancels  out  in  other  places  where  a  mean- velocity  gradient 
or  a  frame  rotation-rate  appears),  r  is  defined  as 

t2=(j)  +(Cr\fi)  '  (48) 

where  Cr  =  2.0.  Away  from  walls,  in  high  Reynolds  number  regions,  this  gives  the  familiar 
definition  r  =  k/e.  Near  the  wall  (or  in  low  Reynolds  number  regions)  the  time  scale  has  a 
lower  bound  in  terms  of  its  viscous  counterpart. 
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The  coefficient  CT  was  chosen  so  that  the  prediction  of  the  algebraic  model  (Section  3.11) 
for  the  energy-normalized  Reynolds  stresses  would  closely  resemble  the  wall  values  given  in 
the  DNS  database  for  channel  flow  (Moser,  Kim,  and  Mansour  1999). 


4.4  Length  Scale 


A  length  scale  is  used  for  the  extension  of  the  k  and  the  u)2  evolution  equations,  with 
the  specific  requirement  that  it  shall  not  vanish  at  a  wall.  For  homogeneous  flows,  a  dimen¬ 
sionally  consistent  form  is 


(4.9) 


where  lh  is  the  homogeneous  length  scale.  Its  wall  value  is  zero.  It  can  be  made  nonzero 
through  the  addition  of  a  Kolmogorov  length  scale,  or  it  can  be  made  nonzero  through  the 
addition  of  an  inviscid  correction.  Both  are  valid  approaches. 

The  selected  form  uses  an  inviscid  correction, 


1  +  Cn 


1  + 


Cn 


\  ,  r  dlh  dlh 
1  +  Cl2d^kd^~k 


(4.10) 


Bu  was  introduced  in  Chapter  3  and  is  representative  of  the  proximity  of  a  wall.  It  is  0(y2), 
so  that  the  form  adopted  renders  l2  finite,  O(y0),  as  a  wall  is  approached. 

The  term  in  the  second  square  brackets  effectively  brings  a  nonlocal  contribution  to  the 
length  scale.  The  idea  is  that  a  given  spatial  location  may  be  affected  by  larger  eddies 
“centered”  in  other  spatial  positions.  In  the  logarithmic  region  of  a  boundary  layer,  for 
instance,  the  length  scale  grows  with  distance  from  the  wall.  It  is  then  possible  to  imagine 
that  the  turbulence  is  locally  affected  by  eddies  larger  than  the  local  eddies.  This  nonlocal 
correction  is  important  for  matching  the  mean  velocity  profile  in  the  logarithmic  region  of 
a  boundary  layer. 


4.5  Turbulent  Kinetic  Energy  Equation 

Substituting  the  triple  decomposition  (Equation  4.1)  into  the  Navier-Stokes  equations, 
Equation  (2.2),  and  applying  successively  the  large-scale  coherent  averaging  (brackets),  and 
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the  time/ensemble  averaging  gives  an  equation  for  mean  velocity  Hi, 


&Ui 

dxi 


=  0, 


dui  _  dui  f_  _  dP*  d2Ui 

Qf  uk  -  —  h  v 


'Ox, 


dxi  dxkdxk 


dRjj 
dxj  5 


(4.11) 

(4.12) 


where  Rij  =  UiUj  +  u"u'-  is  the  total  Reynolds  stress  tensor,  the  third  term  on  the  left  is 
the  Coriolis  acceleration  term,  and  P*  encompasses  the  centripetal  acceleration  term, 

P  —  —  •r(€ijk(ljXk)(eipqQ^>Xq)i  (4*13) 

and  U  is  an  average  conservative  body  force. 

Subtracting  the  equation  for  U{  from  the  equation  obtained  after  the  first  averaging 
process,  ( Ui  -bit*),  gives  the  equation  for  fi*.  Multiplying  by  Ui  gives  the  equation  for  the 
large-scale  kinetic  energy  as  reported  in  RLK02, 


-  H«X»g  - 

'  v  ^  s _  ^  V  j  > 


dt 


DKDT 

mean  Lagrangian  increase 


TKML 

transfer  from  mean 


TKLS 

transfer  to  small  scales 


DKL 

dissipation  by 
large  scales 


d_ 

dx 


spatial  flux 

Similarly  the  small-scale  kinetic  energy  equation  is  (RLK02), 


(4.14) 


-  (~pUj  +  UjUiUi/2  +  Ui{u"u")  +  ~v~^~  j 

'  FKLj  FKS,  ' - - - ' 

FKPj  J  J  FKi/j 


du"u"/2  _  du'fu'f/2  t  „  „w 

. y . +uj . L  . =  +  H«X» 


dxj 


Dk"DT 

mean  Lagrangian  increase 


3  dXn 


duj 
’  dx-j 


Tk"ML 

transfer  from  mean 


—  V 


TKLS 

transfer  from  large  scales 


-p"u"  +  uju"u"/ 2  +  u"u"u"-  + 


Fk"P 


du"  du" 
dxj  dxj 


Dk"L 

dissipation  by 
small  scales 


—V 


du"u"/2 


dx,j 


Fk"i/| 


(4.15) 


spatial  flux 


RLK02  chose  the  total  kinetic  energy  as  the  first  turbulence  scale.  Adding  the  large- 
scale  coherent  kinetic  energy  to  the  small-scale  kinetic  energy,  Equations  (4.14)  and  (4.15) 


respectively,  gives 


dt  +  dxj 


dui  Out 


,  /  n  n\\dUi  du"  du'l 


where  A;  =  w^/2,  w'  =  £  +  it",  and  p'  =  p  +  p". 

The  next  step  can  be  more  controversial.  Clearly  the  term  responsible  for  transfer 
from  the  large  scales  to  the  small  scales,  TKLS,  strictly  cancels  out  and  the  dissipation 
is  represented  by  the  sum  of  the  dissipation  by  large  scales,  DKL,  and  by  small  scales, 
Dk"L.  The  latter  represents  the  bulk  part  of  the  dissipation  in  high  Reynolds  number 
flows.  In  RLK02  the  total  dissipation  is  modeled  as  the  sum  of  DKL  and  TKLS,  with 
the  implied  assumption  that  TKLS  and  Dk"L  balance  each  other.  This  choice  follows  the 
concept  that  what  should  be  modeled  is  the  transfer  of  energy  down  the  cascade,  rather 
than  the  dissipation  at  the  small  scales  directly.  It  has  a  powerful  physical  appeal;  it 
is  attractive  to  model  a  large-scale  process  instead  of  providing  an  ad-hoc  expression  for 
small-scale  phenomena.  Here  the  approach  of  RLK02  is  followed,  the  term  within  brackets  in 
Equation  (4.16)  is  dropped,  and  the  final  equation  for  the  turbulent  kinetic  energy  becomes 


dk  .  dk 
dt  +  dxj 


(4.17) 


where  e,  the  dissipation  rate  of  turbulent  kinetic  energy,  is 


dui  dui 


S  V  1  *  S  V  * 

TKLS  DKL 

transfer  to  small  scales  dissipation  by 
large  scales 


4.5.1  Term-by-Term  Modeling 

The  first  term  on  the  right-hand-side  of  Equation  (4.17)  represents  energy  transfer  from 
the  mean  flow  to  the  large  scales  of  the  turbulence  (best  known  as  production  of  turbulence 
kinetic  energy)  and  does  not  require  modeling.  The  remaining  terms  are  modeled  such 
as  to  resemble  the  near-wall  asymptotic  behavior  of  the  exact  terms  and  to  revert  to  the 
homogeneous  RLK02  form.  Table  4.1  summarizes  the  modifications  done  to  the  model 
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proposed  by  RLK02,  as  well  as  the  near-wall  asymptotic  behavior  of  the  exact  terms.  In 
the  sequence  the  changes  are  briefly  commented.  The  symbol  fd  is  introduced  to  shorten 
the  notation, 


a  =  fffdS. 


(4.19) 


Term 

exact  term 

i 

o 

Current  model 

RLK02 

TKLS 

0(j/3) 

„ ,  ki,~- 
F,k— 

Feku 

Fe  =  Csi&S)1** 

Fe  =  C£(  3fd) 

DKL 

dui  du{ 
dxj  dxj 

0(1) 

vu2 

triple 

product 

0(y3) 

d  fait  Vkldv,iUi/2\ 

dxj  y<Jkrjk  Fe  dxk  J 

- 

pressure 
~l~  transport 

-£;H) 

0(y2) 

Table  4.1:  Summary  of  terms  in  the  large-scale  turbulent  kinetic  energy  equation,  their  asymptotic 
behavior  as  a  wall  at  y  —  0  is  approached,  the  current  model  form,  and  the  original  model  form  in 
RLK02. 

4.5.2  TKLS:  Transfer  to  Small  Scales 


TKLS  =  =  ft*—,  =  ^(3fd)1/6.  (4.20) 

The  current  model  has  the  proper  asymptotic  behavior  as  the  wall  is  approached,  as 
y/k  =  0(y)  and  l  =  0(1)  as  y  — >  0.  More  importantly,  the  nonlocal  modification  to  l2  (the 
0/2  term  in  Equation  4.10)  effectively  changes  Ce  in  different  parts  of  the  flow.  This  mod¬ 
ification  is  necessary  to  match  velocity  profiles  from  DNS  of  channel  flow.  The  drawback 
of  this  modification  is  that  Fe  had  to  be  modified  from  RLK02,  to  incorporate  the  1/6-th 
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power  in  Ce( 3M)1/6.  Fe  is  used  to  shut  down  the  cascade  of  energy  from  large  scales  to 
small  scales  in  rotating  turbulence. 

The  introduction  of  C12  also  implies  that  a ^2  / <1.5  for  proper  behavior  near  the 
edge  of  a  turbulent /non- turbulent  interface  (see  Section  4.8.1),  as  opposed  to  cr^/ak  <  2 
when  Q2  =  0.  As  shown  in  Section  4.8.2,  matching  the  log-law  when  Fe  =  Ce(3fd)  requires 
a^2 /<7k  =  1-9  when  C12  j=-  0.  If  Fe  =  C^fd)1/6  then  a^ijok  —  1.45  when  C12  ^  0. 

In  RLK02  the  specific  value  of  C£  was  found  to  be  immaterial  for  high  Reynolds  number 
flows,  as  it  could  be  absorbed  by  the  other  model  constants.  In  wall-bounded  flows,  the 
value  of  Ce  becomes  an  important  parameter  as  one  bridges  the  near-wall  layer.  Results 
shown  in  the  present  work  were  obtained  with  a  constant  value  for  all  flows, 

C£  =  0.3.  (4.21) 


4.5.3  DKL:  Dissipation  by  Large  Scales 


Expanding  uoj2  —  v&i&i  in  terms  of  velocity  gradients,  and  using  the  continuity  equation 
for  the  fluctuations,  gives  the  identity, 


_  dui  din  _2  d2UjUj 
~  V  dxj  dxj  UlJ  V  dxidxj ' 


(4.22) 


Note  that  DKL  is  strictly  the  homogeneous  dissipation.  It  differs  from  the  “true”  dissipation 
by  gradients  of  the  Reynolds  stresses.  The  true  dissipation  is  positive  definite  and  represents 
the  work  performed  by  a  stress  against  a  deformation. 

The  second  term  in  the  right,  the  second  derivative  term,  is  actually  very  small  compared 
to  the  sum  of  vu> 2  and  TKLS  and  it  will  be  dropped  for  simplicity.  DKL  is  then  modeled 
as 


DKL  =  i/a>2. 


(4.23) 


4.5.4  FK:  Spatial  Flux 

The  transport  term  consists  of  a  turbulent  contribution  and  a  viscous  contribution.  The 
viscous  flux,  FKvj ,  does  not  require  modeling.  When  integration  down  to  the  wall  is  con¬ 
sidered  this  is  a  dominant  term,  0(1)  as  y  — »  0.  Otherwise  it  is  negligible  in  high  Reynolds 
number  regions.  The  flux  of  velocity  fluctuations  and  the  flux  of  pressure  fluctuations  have 
different  near  wall  asymptotic  behaviors.  The  latter  is  0(y 2)  as  y  — >  0,  while  the  former  is 
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0(y4).  The  whole  turbulence-induced  flux  is  traditionally  modeled  as  a  gradient-diffusion 
process.  Cazalbou  and  Bradshaw  (1993)  conducted  a  priori  analysis  of  turbulent  transport 
in  wall-bounded  flows  and  confirmed  good  results  for  the  representation  of  turbulent  trans¬ 
port  in  terms  of  gradient  diffusion  models.  In  DNS  of  channel  flow  at  several  Reynolds 
numbers  (Moser  et  al  1999)  the  pressure  flux  is  significantly  smaller  than  the  velocity  flux 
away  from  the  wall,  and  has  opposite  sign.  Figure  4.1  shows  this  for  the  ReT  =  590  case. 
The  pressure-velocity  correlation  was  available  from  the  DNS  database  and  the  total  flux 
was  integrated  from  its  gradient  given  in  the  turbulent  kinetic  energy  equation  balance  given 
in  the  database  (the  triple  product  follows  by  subtraction).  The  sum  of  both  turbulence 
fluxes  resembles  the  velocity  flux,  and  it  can  be  successfully  modeled  away  from  the  wall 
(for  the  DNS  channel  flow  database)  as 


1 


k  dk 


-p'u 2  +  2  «  -  0.25i?22~7~* 


(4.24) 


*  ,  .  .  .  .  otk  yfk  l  dk 

-p'u'j  +  u'juX/2---rjk——. 


(4.25) 


9  ~  ~^edy 

Figure  4.2  shows  the  total  kinetic  energy  flux  (using  wall  normalization)  versus  distance 
to  the  wall  for  the  DNS  of  channel  flow  at  the  three  Reynolds  numbers,  compared  to  the 
model  above. 

The  final  form  used  in  this  work  aims  to  reproduce  the  right  hand  side  of  Equation  (4.24), 
while  reproducing  the  near-wall  asymptotic  behavior  of  the  triple-product  term, 

1- 

P1  J  J  ‘  •'  CTfc 

Away  from  the  wall  TKLS  dominates  the  dissipation  rate  e ,  so  that  kl^2l  «  k2/e  and 
Equation  (4.25),  applied  to  the  channel  flow,  falls  back  to  Equation  (4.24).  Thus  OLk/vk 
should  be  close  to  2  x  0.25.  Note  that  the  ratio  of  a &  to  <x*.  is  the  important  quantity,  so  that 
Gk  can  be  chosen  as  1,  without  loss  of  generality.  The  model  is  0{y4)  near  the  wall  which 
agrees  with  the  asymptotic  behavior  of  the  triple  product.  The  pressure-transport  term  is 
0(y 3),  so  the  overall  model  transport  near  the  wall  grows  slower  than  the  real  transport. 
This  is  however  a  second  order  effect.  It  is  important  to  get  viscous  transport  much  larger 
than  turbulent  transport  near  the  wall. 


4.6  Large-Scale  Turbulent  Enstrophy 

The  starting  point  to  build  the  large-scale  turbulent  enstrophy  equation  is  the  instan¬ 
taneous  vorticity  equation,  obtained  by  taking  the  curl  of  the  continuity  equation  and  the 
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Navier-Stokes  equation,  Equations  (2.1)  and  (2.2)  respectively, 


duip 

dt 


+  Uj 


dwp 

dxp 

dup 

dxj 


=  0, 


°jp 


WpSpp  +  2nhxPj+l'dxJxk- 


(4.26) 

(4.27) 


The  large-scale  turbulent  enstrophy  equation  is  obtained  following  the  same  procedure 
followed  to  get  the  large-scale  turbulent  kinetic  energy  equation.  In  short,  after  substituting 
the  triple  decomposition  for  the  vorticity,  Equation  (4.3),  the  application  successive  aver¬ 
aging  operations  permits  the  isolation  of  equations  for  the  mean  vorticity,  u>i,  and  for  the 
large-scale  vorticity,  Taking  the  product  of  the  latter  equation  with  u)*,  and  averaging 
produces  the  desired  result.  The  final  exact  equation  for  the  large-scale  turbulent  enstrophy 
is  reproduced  from  RLK02, 


dui&i/ 2  _  duiUJi/2  —-dLJi  ,  .  n 

“at -  +  u>-ei~‘ 

^  ^  V.  ...  ✓  v  ^  ■  ■  ■✓ 


+  CbiCbj 


duj 

dxj 


(4.28) 


DWDT 

mean  Lagrangian  increase 


TWML 

transfer  from  mean 


TWLS 

transfer  to  small  scales 


PLLM 

production  by 
mean  stretching 


+ 


-  -  dui  _  d%k  /  „  du"  \  dwi 

+u’jUid^+UiY^)~  vWi 


dtij 

dxj 


PLLL 

production  by 
self  stretching 


PLML,  PLSS 
production  by 
cross— scale  stretching 


DWL 

dissipation  by 
large  scales 


d 

dxj 


UjU)i<jJi/2  +  uh{u)'M)  + 

v - v - ,  v - - , 

k  FWLj  FWSi 


-V 


dCbiCbi/i 


dx  i 


FWi/i 


spatial  flux 


4.6.1  Term-by-Term  Modeling 

The  “mean  Lagrangian  increase”  and  the  viscous  diffusive  flux  terms  do  not  entail  any 
modeling.  The  other  terms  are  modeled  such  as  to  resemble  the  near-wall  asymptotic 
behavior  of  the  exact  terms,  and  to  have  them  revert  to  the  homogeneous  form  proposed 
by  RLK02.  Note  that  TWML  and  PLML  vanish  for  homogeneous  turbulence.  Table  4.2 
summarizes  the  modifications  imposed  on  the  model  proposed  by  RLK02,  as  well  as  the 
near-wall  asymptotic  behavior  of  the  exact  terms.  In  the  sequence  the  changes  are  briefly 
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commented.  The  symbol  rdf  is  introduced  to  shorten  the  notation, 


rdf  =  njdftfft. 


(4.29) 


Term 

exact  term 

o 

t 

Current  model 

RLK02 

TWML 

<-**>£ 

0{y 2) 

0 

0 

TWLS 

H“K»ig 

0(y) 

n  ~2 

Cq2TU)Z  — 

PLLM 

.  .  dUi 

0(y) 

-.2  icc 

Uf^Wj 

-2 1  9Ui 

uf“w, 

PLLL 

0(y) 

C(j2p(9rdi)u)2~j- 

Cj)2p  (9rdf  )d)3 

+PLSS 

*«> 

PLML 

_  -  dUi 

0(y) 

0 

0 

DWL 

du ){ d&i 
dxj  dxj 

0(1) 

£j2 

C{j2v 

V4 

C^UT 

FWLjj 

+FWSjj 

—  (ntuJM/2 

+^K<>) 

0(y ) 

9  (akr  Vkldu2/2\ 
dxj  Fe  dxk  J 

- 

Table  4.2:  Summary  of  terms  in  the  large-scale  turbulent  enstrophy  equation,  their  asymptotic 
behavior  as  a  wall  at  y  =  0  is  approached,  the  current  model  form,  and  the  original  model  form  in 
RLK02. 
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4.6.2  TWML:  Transfer  from  Mean  Flow 

This  term  vanishes  in  homogeneous  turbulence,  as  the  velocity  gradients  can  be  at  most 
a  function  of  time  so  that  the  spatial  derivative  of  the  mean  vorticity  must  necessarily  van¬ 
ish.  Substituting  the  near-wall  expansions,  Equations  (4.6a-e)  and  (4.7a-c)  into  TWML’s 
definition,  it  is  found  to  be  0(y2)  near  the  wall,  as  indicated  in  Table  4.2.  The  large-scale 
velocity- vorticity  correlation  can  be  separated  into  a  symmetric  part  and  an  antisymmetric 
part.  The  latter  part  has  an  exact  representation  in  terms  of  gradients  of  Reynolds  stresses. 
Consider  the  identity 


Cijk&i&j  —  (WpWp/2)^/*.  (4.30) 

Writing  the  correlation  in  terms  of  its  symmetric/antisymmetric  decomposition, 

UJjUj  —  X-ij  4"  (4.31) 


and  substituting  into  Equation  (4.30)  gives 


W 


—  2  [(^fc^p),p  (ttptfcp/2)^]. 


(4.32) 


There  is  no  exact  representation  for  the  symmetric  part.  Robinson,  Harris,  and  Has- 
san  (1995)  proposed  an  ad-hoc  gradient-diffusion  model  using  an  isotropic  turbulent  eddy 
viscosity.  Generalizing  that  model  to  account  for  anisotropic  Reynolds  stresses  ultimately 
produces, 


Xij  =  i  (uiiUj  +  OjUi) , 


=  C3  (ru  k 


dQj 

dxt 


+  r  ji  k  t 


an 

dxi 


')• 


(4.33) 

(4.34) 


where  C$  is  a  model  constant,  r  is  a  time  scale,  k  is  the  turbulent  kinetic  energy,  and  Qj 
represents  the  mean  flow  vorticity  vector. 

Alternatively,  Bernard  (1990)  conducted  a  Lagrangian  analysis  and  found 


■y-sr-  ^  _  dQj  _ 

ujUi  -  -PRik—  +  Qui—nh, 

where  P  and  Q  are  Lagrangian  integral  time  scales.  It  then  results 

QdRjj 


(4.35) 


__  _  (  ,  dQj  .  dQj\ 

x,>  =  -plr<‘kdi+r»ke^) 


(4.36) 
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Note  that  Bernard’s  form  has  an  extra  term  with  a  different  tensorial  form,  when  compared 
to  Robinson  et  al .  Attempts  to  model  the  terms  containing  L OjUi  were  conducted  here,  but 
no  definitive  conclusion  could  be  reached.  These  terms  are  important  in  the  buffer  layer  of 
a  wall-bounded  flow.  However  there  is  already  uncertainty  on  the  form  of  the  other  terms 
in  the  adequation  in  the  same  region,  so  these  could  be  unintentionally  overcoming  the 
requirement  for  modeling  the  vorticity- velocity  correlation  terms.  As  no  clear  improvement 
could  be  discerned  including  these  terms  so  far,  terms  involving  CjjUi  were  neglected  in  the 
present  analysis. 

4.6.3  [TWLS  -(PLLL+PLSS)] :  Transfer  to  Small  Scales  and  Production 
by  Self-Stretching  and  Cross-Stretching 

Following  RLK02,  this  combination  is  modeled  such  that  in  2D-2C  homogeneous  turbu¬ 
lence  TWLS  does  not  vanish  and  (PLLL+PLSS)  does  vanish.  The  form  chosen  in  RLK02 
is  not  appropriate  in  near-wall  regions,  as  it  predicts  the  sum  of  these  terms  to  be  0(1)  as 
y  —>  0,  contrary  to  the  exact  term.  It  is  important  to  fix  this  anomaly  as  the  near- wall  bal¬ 
ance  should  be  restricted  only  to  the  viscous  diffusion  term  and  the  large-scale  dissipation 
term.  To  reconcile  the  RLK02  form  with  the  desired  near-wall  asymptotic  behavior  it  is 
altered  to 

(«7) 

RLK02  found  that  when  Fe  =  (3fd),  C^TjCt  =  25/14  and  C^P/Ce  =  20/21  matched 
decay  rates  of  homogeneous  turbulence  in  fixed  and  rotating  frames  of  reference.  The  anal¬ 
ysis  of  RLK02  was  revisited  to  match  the  same  decay  rates  with  the  updated  Fe  =  (3fd)1//6 
(for  a  small  wavenumber  part  of  the  spectrum  proportional  to  the  square  of  the  wavenumber 
as  in  the  selection  of  RLK02  chosen),  resulting  in 

=  2.56,  =  1.73.  (4.38) 

C£  ce 

4.6.4  PLLM:  Production  by  Mean  Stretching 

The  model  remains  the  same  used  in  RLK02.  Its  essence  is  the  fact  that  the  normalized 
correlation  Ui&j/u)2  is  intimately  related  to  the  normalized  circulicity  tensor,  as  can  be 
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inferred  from  Chapter  2, 


^r-^-dUi 


-2  ccc  dUi 


(4.39) 


4.6.5  PLML:  Production  by  Cross- Stretching 


This  term  is  identically  zero  in  homogeneous  turbulence.  This  result  follows  from 


.  dui  _  ( diii 
(Joe  j  y  c/  Xj 


duj\ 

dxi) 


_  dilj 
+  u>id. V 


—  Ckj&k&i  + 


duiiUj 

dxi 


duiiUj 
dxi  ’ 


(4.40) 

(4.41) 


which  vanishes  in  homogeneous  situations.  To  get  the  last  expression  note  that  the  large- 
scale  vorticity  is  divergence-free.  Closure  for  this  term  again  requires  a  model  for  OiUj.  As 
in  the  TWML  term,  this  will  also  be  kept  zero  for  the  current  model. 


4.6.6  DWL:  Dissipation  by  Large  Scales 


The  choice  used  in  RLK02  is  not  appropriate  in  near-wall  regions.  It  is  made  finite  at 
the  wall  by  introducing  the  length  scale  in  its  definition, 


duj  dujj 
dxj  d  Xj 


n 

:  Ctfv  V  p  ■ 


(4.42) 


RLK02  find  C —  5/6  to  match  the  turbulent  kinetic  energy  decay  rate  at  small  turbulent 
Reynolds  numbers  (final  period  of  decay). 


4.7  Wall  Boundary-Conditions 

The  boundary  conditions  used  at  a  solid  wall  refer  only  to  the  turbulent  kinetic  energy, 
k.  It  is  constrained  to  satisfy  its  near-wall  asymptotic  behavior.  There  is  no  boundary 
condition  that  needs  to  be  satisfied  by  the  large-scale  turbulent  enstrophy.  It  adjusts  itself 
near  the  wall  to  ensure  that  the  proper  behavior  of  the  turbulent  kinetic  energy  is  enforced. 

The  boundary  conditions  are  then  applied  to  the  couple  system  k-oj 2  as 

Qk 

*  =  0,  ^=0,  (4.43) 


where  xn  is  the  direction  normal  to  the  wall. 
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4.8  Compatibility  with  Canonical  Flow  Regions 

Enforcing  simultaneous  compatibility  with  a  turbulent /non- turbulent  interface  and  a 
logarithmic  region  determines  the  ratio  of  the  numerical  coefficient  of  the  transport  term 
for  the  scale  equations.  It  also  reinforces  model  behavior  consistent  with  observed  log-law 
behavior,  as  well  as  lack  of  spurious  solutions  near  the  edge  of  a  turbulent/non-turbulent 
interface. 


4,8.1  The  Edge  of  a  Turbulent  Region:  Turbulent /Non- Turbulent  Inter¬ 
face 

The  edge  of  a  turbulent  region  should  not  be  important  as  far  as  increased  mixing  due 
to  turbulence  is  concerned,  as  very  little  happens  there.  It  is  important  as  it  may  relate  to 
the  entrainment  of  non-turbulent  flow  into  a  turbulent  flow  region.  In  modeling  it  relates 
to  the  sensitivity  of  the  solution  inside  a  turbulent  layer  to  the  conditions  outside.  The 
objective  here  is  to  find  potential  constraints  regarding  model  constants,  in  order  to  avoid 
spurious  solutions. 

A  turbulent  front  is  investigated  in  an  effort  to  ensure  the  model  is  consistent  with  a 
well-behaved  solution  in  this  region.  Cazalbou,  Spalart,  and  Bradshaw  (1994)  and  Cazal- 
bou  and  Chassaing  (2002)  discuss  the  essential  physics  associated  with  the  presence  of  a 
front.  In  short,  a  sharp  interface  between  turbulent  and  non-turbulent  fluid  is  expected  so 
that  turbulence  statistics  should  vanish  at  a  finite  distance  instead  of  extending  into  the 
shear-free  region.  Additionally  they  notice  that  budgets  of  modeled  and  real  turbulence 
kinetic  energy  should  reduce  to  a  balance  between  turbulent  diffusion  and  mean  convective 
transport.  The  same  applies  to  the  modeled  second  scale  equation.  A  simple  model  problem 
to  investigate  the  front  region  is  proposed;  the  unsteady  nonlinear  diffusion  problem.  In 
one  space  dimension,  particularized  to  the  current  model, 

dU  _  dRu 
dt  dy  ’ 
dk  _  d  fvr  dk\ 
dt  ~  dy  V<7fc  dy)  ’ 
du)2/2  _  d  /  vt  d£j2/2\ 
dt  ~  dy  V<7£2  dy  )  ’ 


(4.44) 

(4.45) 

(4.46) 
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where 


vT  =  2C„  r22-^-.  (4.47) 

Near  the  front  the  length  scale  and  the  dissipation  of  the  turbulent  kinetic  energy,  Equations 
(4.10)  and  (4.18),  become 


l2  = 


u4. 


s  —  Fek 


1  +  C12 


d(Vk/ Cj)  d(Vk/to) 


dy 


dy 


(4.48) 

(4.49) 


Cazalbou,  Spalart,  and  Bradshaw  (1994)  find  solutions  in  terms  of  the  distance  to  the 
edge, 


z  =  ct  -y,  (4.50) 

where  c  is  the  propagation  velocity  of  the  front.  Following  Paullay  et  al.  (1985),  the  solution 
near  the  edge  is  found  in  terms  of  a  power  expansion  in  terms  of  the  distance  to  the  edge, 

U  =  Auzn ,  k  =  Akzp,  Cj  —  Awzq,  (4.51) 

where  n,  p,  g,  Au,  A^)  and  A&  should  all  be  positive.  The  model  problem  and  its  solution 
assume  that: 

•  the  turbulence  generated  in  the  core  of  the  shear  layer  spreads  with  a  finite  velocity 
into  the  non-turbulent  fluid  —  thus  displaying  a  sharp  interface  at  a  finite  distance; 

•  the  fronts  of  the  mean  velocity,  turbulence  kinetic  energy,  and  the  second  scale  always 
coincide; 

•  the  budget  of  each  scale  equation  reduces  to  a  balance  between  the  time  rate  of  change 
and  turbulent  diffusion  in  the  vicinity  of  the  front. 

Substituting  Equations  (4.51)  into  Equation  (4.45)  (turbulent  kinetic  energy)  gives 
Atcpzv~l  =  ~  Wt±tl  ~  (4.52) 


Comparing  powers  of  z, 


Q  =  Ip  —  I- 


(4.53) 
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Comparing  coefficients, 

aI/2VQ2Cu^(j>  -  4)2r22p  (  ) 

2AlFeak  '  1  J 

Substituting  Equations  (4.51)  into  Equation  (4.46)  (large-scale  enstrophy),  using  Equa¬ 
tions  (4.53)  and  (4.54),  and  comparing  coefficients  gives  p, 

4  (?k 


P  = 


3 <Tk  ~ 


Thus 


2  Q  = 


40^2 


(4.55) 


(4.56) 


3ak  ~  2^2 

For  compliance  with  the  current  analysis  Equations  (4.55)  and  (4.56)  require  a^/ok  >  3/2. 

In  order  to  analyze  the  [/-momentum  equation,  the  functional  form  of  #12  has  to  be 
specified.  Assuming  ( dU/dy)r  approaches  zero  at  the  front,  the  current  ASBM  (Chapter  3) 
produces  the  limiting  form  of  R12  (following  an  analysis  similar  to  Haire  2003), 


_  1  ,dUk 

&12  =  - . 

3ao  ay  e 

Near  the  front  the  dissipation  rate,  e,  and  the  time  scale,  r  =  fc/s,  become 

A^TkAlFez-^ 


e  = 


VW(P~4)2 


k  =  p 

T  e  4AIF£ 


(4.57) 


(4.58) 


(4.59) 


Substituting  Equations  (4.51)  into  Equation  (4.44)  (mean  momentum),  using  Equations 
(4.53)  and  (4.54),  and  comparing  coefficients  gives  n,  the  exponent  of  the  velocity  expansion, 

n  =  eapCm  (4  60) 


°k 


n  will  be  positive  for  p  >  0. 

The  condition  to  be  satisfied  in  order  for  ( dU/dy)r  — ►  0  can  be  found  by  expressing  it 
in  terms  of  the  expansions  above, 


dU  y/^kAuVUi2^/(p  -  4 )2nzn  p 

dyT~~  4AlFe 


It  will  approach  zero  if  (n  —  p)  >  0, 

n —p =p ^ 


6apCur22  —  &k 

°k 


)■ 


When  ( dU/dy)r  — *  0,  it  follows  that  7*22  —►  1/3,  so  the  constraint  becomes 

2g&oCi/  ^  o h . 


(4.61) 


(4.62) 


(4.63) 


Compatibility  Requirement  Summary 


Summarizing,  the  near-edge  analysis  yields  the  following  constraints,  which  should  be 
satisfied  in  order  to  achieve  a  well-behaved  solution: 

•  CT^/cjk  <  3/2  (Equations  4.55,  4.56,  and  4.60); 

•  2aoCu  >  <7k  (Equation  4.63). 

The  values  used  in  the  present  work,  do  =  1.6,  Cv  —  0.25,  and  =  1  violate  the 
constraint  above.  As  a  compromise  solution,  ao  could  itself  be  made  a  function  of  the 
normalized  strain  rate.  It  would  be  made  higher  at  zero  applied  strain  rate,  and  then  it 
would  taper  off  towards  the  optimized  value  above.  For  example, 

ao  =  1.6  +  a'0  exp(-ao'52r2).  (4.64) 

This  work  is  concentrated  on  wall-bounded  flows  which  do  not  have  free-stream  boundaries. 
It  is  then  futile  to  optimize  the  amplitude  (a[,)  and  delay  (a^)  of  the  above  modification 
here.  Any  optimization  should  be  done  with  a  flow  where  it  matters.  For  instance,  Haire 
(2003)  found  that  the  spreading  rates  of  self-similar  free  shear  flows  varied  inversely  with 
ao-  The  optimization  is  left  as  a  recommendation  for  future  work. 


4.8.2  The  Log  Layer 


The  log  layer  is  investigated  in  detail  to  ensure  that  the  model  is  consistent  with  observed 
behavior  of  the  mean  flow  and  the  turbulence  quantities.  In  the  log  layer  of  a  high-Reynolds- 
number  flow,  the  scale  equations  become 


dk  n  .  _  .du  „uVk  d 
di~  0 -<-*>’>«*- +  5; 

^2/2_n_,-.2  tccW  nrfiVk.d 


dt 


=o=^m--c^  +  - 


=  ^ 


1  +  C12 


(dVk/u,\‘ 

l  dy  ) 


(a/t  y/k  l  dk  \ 

)  ’ 

(  ak  Vkldd)2/ 2\ 

(a*/22  Fe  dy  J  ’ 


(4.65) 

(4.66) 

(4.67) 


where  Co  =  \cC)2T  -  (9ryd^/jf)Css.pJ  and  ak  =  1  (Section  4.5.4). 
In  a  high-Reynolds-number  log  layer 


dU  __  uT 
dy  Ky ’ 


-Rn 

k 


ut 


k 


—  2(-ri2)  =  ai, 


k  —  constant. 


(4.68) 
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where  u*  is  the  wall  shear-stress.  cq  =  0.3  is  a  well-accepted  value  (showing  up  a sCM  =  0.09 
in  most  eddy-viscosity  models).  The  Von  Karman  constant,  «,  is  somewhat  more  contro¬ 
versial.  Traditionally  k  =  0.41  has  been  used.  However  there  are  recent  indications  that 
high  Reynolds  numbers  boundary  layers  may  be  better  represented  by  smaller  values,  say 
k  =  0.38  (Osterlund  et  al  2000).  The  normalized  turbulence  structure  is  also  constant 
in  the  log  layer.  Kim  (1992)  gives  the  normalized  structure  tensors  in  a  channel  flow  at 
Rer  =  385.  Even  though  it  is  a  fairly  low  Reynolds  number,  and  k  varies  throughout  the 
layer  (as  does  Rn),  the  normalized  structure  is  fairly  constant  with  representative  values 
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0 
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and  it  follows  that 


9  ri:jd?kfg  =  0.65,  SUjdji  =  0.90.  (4.70) 


In  a  high  Reynolds  number  log  layer,  the  Reynolds  shear  stress  is  roughly  constant, 
(—R12)  =  w?.  By  virtue  of  constant  turbulent  kinetic  energy,  the  transport  term  in  its 
equation  vanishes,  resulting  in  a  balance  between  production  and  dissipation,  P  —  e.  This 
implies  the  following  log-law  relationships, 


Z  = 


Vk 

l 

FEk 3/2 


2 

4— 


ut 


ut 


-ny  = 


Ky  ny 
F€k 

(— 2n2)3/2 


v  =  ky- 


(4.71) 

(4.72) 


Substituting  Equation  (4.72)  into  Equation  (4.67)  and  solving  for  ui. 


Now  substituting  Equations  (4.72),  (4.73),  and  (4.68)  into  Equation  (4.66), 


(4.73) 


y 


Ky 


4m- 

nit  nl  Un  oq2  Fe  dy  \  yA  J 


dy 


0, 


+«."!»:  2k  =  0. 

k  F£k  (7^2  |2riar 


fee 

J 12 


thus 


l  _  2  (n2): 

(Tq2  K2ak  1'22 


(4.74) 
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Enforcing  compatibility  with  the  log  law  constrains  the  value  of  Conversely,  esti¬ 
mating  o^i  according  to  Equation  (4.74)  dictates  the  predicted  value  of  k  for  the  model 
Note  that  both  Ce  and  Cn  are  absent  from  the  constraint,  indicating  that  the  slope  of  the 
log  law  in  high  Reynolds  number  flows  is  independent  of  them.  However,  both  constants 
affect  the  solution  in  the  near-wall  region,  thus  affecting  the  constant  (intercept)  in  the 
logarithmic  law  for  the  mean  velocity. 

Compatibility  Requirement 

Using  k  =  0.38  and  the  structure  given  by  Equations  (4.69)  and  (4.70),  with  a k  =  0.25 
(crk  =  1),  then  if  Fe  =  C£(3fd)  as  in  RLK02  (with  corresponding  Co),  results  in  cr^  =  1.93. 
This  violates  the  edge  analysis  requirement  (7^2 /ak  <  3/2.  For  compliance  with  both  the 
log  law  and  the  front  analysis,  Fe  was  modified  in  the  current  work  (and  corresponding  Co). 
Choosing  it  to  be  Fe  =  Ce(3fd)1/6,  results  in 

<7-2  =  1.45.  (4.75) 

bringing  <7^2  into  compliance  with  both  the  log  law  and  the  turbulence  front  constraints. 


4.9  Summary  of  the  Scale  Equations 


To  set  up  model  scale  equations  with  a  firm  background,  standard  Reynolds  decomposi¬ 
tion  (mean  and  fluctuating  parts)  was  abandoned  in  favor  of  a  triple  decomposition  (mean, 
large-scale  fluctuations,  and  small-scale  fluctuations).  The  latter  distinguishes  between  the 
large  turbulence  scales,  which  are  responsible  for  large-scale  mixing  and  which  ultimately 
set  the  dissipation  rate  of  the  turbulent  kinetic  energy,  and  the  small  (dissipative)  scales, 
which  mostly  react  to  the  energy  sent  down  to  the  small  scales  from  the  large  scales. 

The  evolution  equation  for  the  turbulent  kinetic  energy  is  not  much  different  from  stan¬ 
dard  forms  used  in  the  literature.  The  greatest  difference  is  the  estimation  of  the  dissipation 
in  terms  of  the  large  scale  enstrophy  of  the  turbulence.  The  evolution  equation  is 


dk  _  dk  __  dui  d 

dt  Uj  dxj  J  dxj  €  +  dxj 


ak  Vkl 
v5jk  +  ^rjk  ^ 


(4.76) 


where  e,  the  dissipation  rate  of  turbulent  kinetic  energy,  is 

e  =  Fek^  +  vw2,  Fe  =  Ce{  3fd),  M  = 


(4.77) 
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The  large-scale  turbulent  enstrophy  equation  is  the  recent  novelty  that  incorporates  a  more 
fundamental  basis.  It  is  given  by 


du>2/2  ,  dw2/2  _  r.2fCCdUi  r„  _  ,  ,2Vk 


dt 


4 -Uj 


dxj 


~  ~  [Cu2T  -  (9rdf)C'(Z,2p]  Lbl—r - CQ2U  v 


+ 


±( 
dxj  y 


a*  Vkl 
v5ik  +  —jjk-jr 


du2/2\ 


where  rdf  =  njdfjft  • 

The  length  scale  is  evaluated  from 


i2  = 


il 


1  4-  Cn 


14- 


Cn 


(1  -Bkky 


1  ,  ri  ®lh 
Cl2dxk  dxk 


uP 

P 


V  —  _ 
lh  ~~  ^9- 


(4.78) 


(4.10) 


The  model  constants  for  the  scale  equations  were  optimized  as  in  RLK02,  and  for  match¬ 
ing  the  channel  DNS  database  of  Moser,  Kim,  and  Mansour  (1999). 


Ce  —  0.3,  C^j'jCe  —  2.56,  C(j2p/Ce  —  1.73,  C^2V  —  5/6, 
Cn  =  0.1,  C12  =  20,  ak  =  0.25,  ak  =  1,  0^2  =  1.45. 


(4.79) 


The  algebraic  equations  developed  in  Chapter  3  respond  to  the  mean  deformation  nor¬ 
malized  with  a  turbulence  time  scale,  t,  defined  as 


4‘)2+M2' 


(4.8) 


where  CT  =  2.0  (strictly  r  is  needed  only  for  normalizing  the  mean  strain  rate  in  Equa¬ 
tion  3.15  —  it  cancels  out  in  other  places  where  a  mean- velocity  gradient  or  a  frame  rotation- 
rate  appears). 

At  a  solid  wall,  the  scale  equations  are  subjected  to  the  boundary  conditions 

k  =  0,  J|  =  0,  (4.80) 

where  xn  is  the  direction  normal  to  the  wall. 
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Figure  4.1:  Turbulent  kinetic  energy  fluxes  in  fully-developed  channel-flow  at  Rer  =  590  (DNS 
by  Moser  et  al.  1999).  Dotted  line:  (1  /p)j/uf2.  Dashed  line:  u,2u'iu,ilel.  Solid  line:  sum, 
(l/p)p'u'2  +  u,2v!iu,il2.  Dot-dashed  line:  model,  -0.2bR22{k/e)(dk/dy). 


Figure  4.2:  Comparison  of  DNS  and  modeled  turbulent  kinetic  energy  fluxes  in  fully-developed 
channel-flow  at  Rer  =  180,  395,  and  590.  Solid  line:  DNS  by  Moser  et  al.  (1999).  Line:  model 
approximation,  ctk/vk  =  0.25. 


Chapter  5 


Model  Results 


Chapter  3  explored  the  response  of  the  algebraic  structure  equations  (uncoupled  from  the 
turbulence)  to  different  homogeneous  normalized  deformation  rates.  This  chapter  explores 
the  complete  ASBM  with  emphasis  on  inhomogeneous  flow  conditions.  The  algebraic  de¬ 
scription  of  the  turbulence  stresses  and  structure,  given  in  Chapter  3,  is  now  combined  with 
the  scale  equations  introduced  in  Chapter  4  to  form  a  complete  set  of  equations  amenable 
to  numerical  integration. 

Initially  the  equilibrium  response  of  the  full  set  of  model  equations  to  homogeneous  shear 
flow  in  a  rotating  frame  is  analyzed.  Pettersson  Reif,  Durbin,  and  Ooi  (1999)  explored  this 
flow  using  it  as  a  guide  to  calibrate  model  constants  and  to  ensure  a  physically  realistic 
model.  The  behavior  of  the  ASBM  is  assessed  for  the  same  flow  to  help  further  guide  the 
choice  of  model  constants. 

Later  the  ASBM  is  assessed/optimized  in  a  set  of  inhomogeneous  flows.  These  corre¬ 
spond  to  a  family  of  pressure-driven  fully-developed  channel  flows,  in  presence  of  distinct 
directions  of  frame  rotation;  spanwise  frame  rotation  aligned  with  the  mean  flow  rotation, 
and  streamwise  frame  rotation  aligned  with  the  principal  flow  direction  (that  of  the  pressure 
gradient). 

The  channel  flow  may  seem  deceptively  simple.  However  it  is  a  stringent  test  to  evaluate 
the  model  capability  to  cope  with  the  integration  of  the  complete  model  equations  through 
the  buffer  layer  to  the  wall.  It  is  already  a  challenge  to  propose  a  physically  acceptable 
and  realizable  turbulence  model,  adequate  for  high-Reynolds-number  flows.  Proposing  a 
realizable  structure  model  adequate  for  integration  through  wall-affected  low-Reynolds- 
number  regions  proved  to  be  quite  difficult. 

Furthermore,  the  presence  of  the  distinct  modes  of  frame  rotation  allows  the  investi¬ 
gation  of  the  model’s  response  to  rotation  in  perhaps  the  most  critical  fashion,  when  all 
secondary  effects  are  solely  due  to  the  turbulence.  In  more  complex  flows,  other  terms 
representative  of  the  mean  motion  (convective,  Coriolis,  body  force,  and  especially  pressure 
gradient  terms)  can  also  be  active  in  the  momentum  equations,  sharing  some  of  the  job  of 
representing  distinct  physical  phenomena.  Then  inaccuracies  in  the  turbulence  model  itself 
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may  be  overshadowed  due  to  smaller  importance  of  the  Reynolds  stress  in  the  momentum 
equations. 


5.1  Homogeneous  Shear  in  a  Rotating  Frame 


Homogeneous  turbulence  subjected  to  shear  in  a  spanwise-rotating  frame  continuously 
evolves  in  time.  Nevertheless  there  exists  an  equilibrium  limiting  state  for  appropriately 
normalized  variables  such  as  the  energy-normalized  Reynolds  stress  and  structure-dimen¬ 
sionality  tensors,  as  well  as  the  ratio  of  the  mean  to  turbulence  time  scales  and  the  ratio 
of  production  to  dissipation  of  the  turbulent  kinetic  energy  (see  for  example  Speziale  and 
Mac  Giolla  Mhuiris  1989). 

The  turbulence  grows  or  decays  depending  on  the  value  of  the  parameter  77/  =  — 20{/5, 
the  ratio  of  frame  vorticity  to  mean  flow  vorticity.  When  the  mean  flow  rotation  is  in  the 
same  direction  as  the  frame  rotation,  77/  <  0,  the  turbulence  experiences  decay  whenever  the 
magnitude  of  77/  is  larger  than  some  small  threshold.  In  the  counter-rotating  case,  when  the 
mean  flow  rotation  is  opposite  to  the  mean  flow  rotation,  the  turbulence  first  experiences 
more  vigorous  growth  in  the  range  0  <  77/  <  1,  and  then  again  it  experiences  decay  when 
77/  grows  somewhat  larger  than  1.  This  coincides  with  the  behavior  of  the  2D  RDT  solution 
for  the  same  flow  given  in  Chapter  3. 

Second-moment  closure  models  display  a  bifurcation  between  two  stable  equilibrium 
states.  The  bifurcation  corresponds  to  a  change  between  exponential  and  algebraic  time 
dependent  growth  of  the  turbulence.  Pettersson  Reif,  Durbin,  and  Ooi  (1999)  introduced 
the  use  of  this  bifurcation  diagram  for  the  improvement  of  turbulence  models,  both  as  a 
guide  to  calibrate  model  constants  and  to  ensure  a  physically  realistic  model. 

The  equilibrium  response  of  the  full  set  of  model  equations  from  the  ASBM  will  be 
investigated  now.  The  mean  flow  is  given  by  Ui  =  {t/ (?/),  0, 0}  in  a  coordinate  system 
Xi  —  {#,  2/,  2},  and  the  frame  rotation  rate  vector  is  given  by  Q{  =  {0, 0,  D3  }.  In  the  following 
the  notation  S  =  dU/dy  is  used  often. 

As  the  turbulence  is  homogeneous,  gradients  of  the  Reynolds  stresses  vanish,  and  the 
mean  flow  remains  constant,  while  the  turbulence  scales  evolve  in  time.  For  the  scales  used 
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with  the  ASBM  the  time  evolution  equations  are 


dk  n  dU 

—  —  — i^12— - e, 

dt  dy 

du;2/2  o.  dU  .  t  .  ,  _  lw2\^ 

— =  w2/i2-^-  -  [Q,2T  -  (9rjfc/fcjdjf)C'i,2p]  w2— , 


(5.1) 

(5.2) 


where  £  and  Z,  respectively  the  dissipation  rate  of  turbulent  kinetic  energy  and  the  length 
scale,  are 


£  = 


Ce(3fijdji)k 


s/k 

l  ’ 


and  the  model  constants  are 


(5.3) 


Ce  =  0.3,  Cspt/Ce  =  2.56,  C&P/Ce  =  1.73. 


(5.4) 


The  time  scale  used  to  normalize  the  mean  deformation  and  provide  input  to  the  ASBM  is 
defined  as 


The  Reynolds  stresses  and  the  turbulence  structure  are  then  computed  with  the  algebraic 
equations  summarized  in  Section  3.11. 

Pettersson  Reif  et  al  looked  at  two  different  measures  of  model  response.  A  bifurcation 
of  the  time  scale  ratio  e/(Sk)  in  77/  space,  and  the  behavior  of  the  ratio  of  production  to 
dissipation  of  the  turbulent  kinetic  energy.  The  existence  of  the  bifurcation  point  depends  on 
a  specific  model  predicting  7*12  oc  O^Sk/e)-1)  as  Sk/e  —>  00.  The  ASBM  however  predicts 
7*12  oc  0((Sk/e)~1/2)  so  it  does  not  predict  a  bifurcation.  This  is  not  necessarily  bad,  as 
the  bifurcation  is  really  a  feature  of  turbulence  models;  there  is  no  evidence  to  what  should 
be  the  behavior  of  real  turbulence.  Thus  the  most  critical  parameter  should  be  the  ratio 
of  production  to  dissipation  of  the  turbulent  kinetic  energy,  P/e,  as  that  dictates  if  the 
turbulence  is  able  to  sustain  itself  or  not. 

The  ASBM  has  this  turbulence  growth  or  decay  dependence  on  77/  built  into  it,  as  the 
structure  parameters  were  set  in  Chapter  3  to  mimic  RDT  behavior  (when  applicable). 
RDT  predicts  turbulence  growth  when  0  <  77/  <  1,  and  turbulence  decay  otherwise. 

Figure  5.1  shows  P/e  for  the  complete  ASBM  model.  If  P/e  >  1  there  is  turbulence 
growth,  otherwise  if  P/e  <  1  the  turbulence  decays.  The  ASBM  prediction  for  the  threshold 
between  growth  and  decay,  rjf  =  {-0.27, 1.96}  happens  at  values  slightly  further  away  from 
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the  region  0  <  77/  <  1  than  suggested  by  Pettersson  Reif  et  al ,  who  quote  the  quasi-linear 
SSG  model  (Speziale,  Sarkar,  and  Gatski  1991)  to  have  P/e  =  1  at  rjf  =  {—0.19, 1.08}.  If 
the  ASBM  constants  were  chosen  so  as  to  have  P/e  =  1  at  77/  closer  to  1,  the  slope  of 
the  mean  velocity  for  the  case  ( Rer  =  180,  Ro  =  0.77)  would  have  improved  somewhat  (see 
Section  5.4),  however  the  magnitude  of  the  turbulent  intensities  for  all  rotating  cases  would 
degrade  significantly. 


5.2  Channel  Flow  Equations 


In  the  sequence  the  mean-flow  and  turbulence-evolution  equations  are  particularized  for 
the  channel  flow  in  a  (possibly)  rotating  frame.  The  equations  are  then  integrated  and 
the  resulting  variables  are  plotted  for  a  plain  fixed  frame  channel  flow,  a  spanwise-rotating 
channel  flow,  and  a  streamwise-rotating  channel  flow. 

The  mean  flow  is  given  by  U{  —  {U(y),0>W(y)}  in  a  coordinate  system  X{  =  {x,y,  z} 
where  y  is  the  wall  normal  direction  (the  sole  direction  of  inhomogeneity),  x  is  the  streamwise 
direction  with  associated  velocity  U,  and  z  is  the  spanwise  direction  with  associated  velocity 
W.  The  frame  rotation  rate  vector  is  given  by  Of  =  {ft{,  0,  Q{}.  The  wall-normal  mean 
velocity  necessarily  vanishes  by  continuity  for  a  fully  developed  channel  flow  with  zero 
velocity  at  the  walls. 

The  mean  momentum  equation  in  the  streamwise  direction  (P*  as  given  by  Equation  4.13 
encompasses  the  centripetal  acceleration), 

dU__dJ»  d2U  dRn 
dt  dx  +  ^  dy2  dy 

can  be  integrated  from  the  wall  at  y  =  0  to  the  wall  at  y  —  2 h,  where  h  is  the  half  channel 
height,  to  yield  (P12  vanishes  at  both  walls), 


dU 

-^(2k) = - 


+  V 


y— 0 


dU 

dy 


y-2h 


—  ur\y=^  ^r\y=^h  —  2ur, 

dP*  u2 
dx  h  ’ 


+  — 
y=o  P 


y-2h 


where  rw  represents  a  wall  shear  stress,  and  the  friction  velocity,  ur  —  yjrw/p ,  may  be  dif¬ 
ferent  on  different  walls,  as  is  the  case  in  the  spanwise-rotating  channel.  In  wall  coordinates, 

1  _  1 

dx+  Rer  ’ 


(5.6) 
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where  P+  =  P* /u2,  x+  =  uTx/v,  and  ReT  =  h+  =  uTh/v. 

The  complete  equations  for  the  mean  flow  and  for  the  turbulence  are 

dU  dP*  d2U  dRn 

dt  dx  + V  dy 2  dy  ’ 

dW  n  d?W  dR23 

dt  V  dy 2  dy  ' 

dk  n  „  dU  „  dW  d  (\  ak  Vkl1dk\ 
dt~°-  Rl2dy  Rn  dy  £  +  dy[[1' +  oj22  Fe  \dyj' 

^2  =  0  =  +  /g2£)  -  [C^T  -  flSHXW]  ^  -  C<to  v\ 

d  (\  ak  Vkl]du2/2\ 


+  37  2 


where  rdf  =  and  €  and  Z,  respectively  the  dissipation  rate  of  turbulent  kinetic 

energy  and  the  length  scale,  are 


6  =  Ce(  fd  =  fijdjl 


2  _  (fc/62)  r  Cn 
1  +  Cn[  (1  -Bii)2 


1  +  Cl2 


dy/k]tf_ 

dy 


and  the  blockage  is  given  by 


Bij  =  6^5 j2$, 

,cP$  . 


The  time  scale  used  to  normalize  the  mean  deformation  and  provide  input  to  the  ASBM 


is  defined  as 


= (!) +  (°V!)  • 


with  CT  =  2.0.  The  Reynolds  stresses  and  the  turbulence  structure  are  then  computed  with 
the  algebraic  equations  summarized  in  Section  3.11. 

The  remaining  model  constants  are  posted  here  again  for  completeness, 


Ce  =  0.3,  C&2T  =  2.56,  Cq2P  =  1.73,  C&v  =  5/6,  (4.79) 

Cn  =  0.1,  Cn  =  20,  ak  =  0.25,  ak  =  1,  <7*2  =  1.45, 
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and  Lurjv  —  23. 

The  model  equations  are  subjected  to  the  wall  boundary-conditions 


17  =  0,  k  =  0, 


#  =  1. 


(5.13) 


Combining  the  scale  equations  with  the  algebraic  equations  for  the  turbulence  structure 
given  in  Section  3.11  gives  the  complete  set  of  equations  solved  to  obtain  the  results  pre¬ 
sented  in  this  chapter.  The  solutions  depend  on  two  parameters  that  completely  characterize 
these  fully  developed  channel  flows;  the  friction  Reynolds  number,  Rer, 


R&J-  — 


uTh 

v 


(5.14) 


and  a  rotation  number,  defined  differently  for  the  spanwise-rotating  cases  and  for  the 
streamwise-rotating  cases. 


5.3  Fixed  Frame  Channel  Flow 

Here  fi{  =  =  0.  Figure  5.2  shows  mean  velocity  profiles  in  wall  coordinates  (normal- 

ized  by  the  wall  shear  stress  and  viscosity)  for  a  series  of  friction  Reynolds  numbers  in  a 
fixed  frame  channel  flow.  Comparisons  are  made  with  the  DNS  of  Moser,  Kim,  and  Mansour 
(1999)  at  Rer  —  {180,395,590}.  Two  distinct  log  laws  are  also  shown.  The  traditional  one 
given  by  Coles  and  Hirst  (1969),  and  the  revised  one  given  by  Osterlund,  Johansson,  Nagib, 
and  Hites  (2000).  The  ASBM  constants  Ce,  C*2,  and  L  were  optimized  for  matching  the 
DNS  at  Rer  =  {395,590}  ( C€  are  C12  not  really  independent  —  given  one  there  is  a  single 
value  for  the  other  that  optimizes  the  mean  velocity  profile).  The  end  result  is  that  the 
DNS  at  Rer  =  180  is  somewhat  underpredicted  while  the  log  law  is  slightly  over  predicted. 

Figure  5.3  shows  the  turbulence  intensities  for  the  fixed  frame  channel  at  Rer  =  590. 
The  peak  kinetic  energy  near  the  wall,  which  can  be  inferred  from  the  intensities,  is  also 
optimized  by  the  selection  of  the  model  constants  Ce,  C/2,  and  L.  The  anisotropy  predicted 
in  the  log  region  is  a  testament  to  the  accuracy  of  the  ASBM  in  this  case.  As  the  centerline  is 
approached,  a  degradation  of  the  results  is  seen:  while  the  DNS  still  predicts  some  anisotropy 
there,  the  ASBM  predicts  isotropy  as  a  result  of  the  zero  deformation  rate  at  the  centerline. 
This  is  a  flaw  with  the  equilibrium  approach.  Near  the  centerline  the  diffusion  of  these 
quantities  becomes  important  and  the  local  equilibrium  assumption  embedded  in  the  model 
no  longer  applies.  In  fact  the  balance  in  the  evolution  equations  for  structure  quantities 
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( e.g .  Reynolds  stresses  or  structure  dimensionalities)  near  the  centerline  lies  in  dissipation 
and  turbulent  transport.  However  this  region  should  not  play  a  major  role  in  the  dynamics 
of  the  complete  flow. 

Figure  5.4  shows  assorted  model  variables  for  the  Rer  —  395  case. 

Figure  5.4(a)  shows  the  normalized  Reynolds  stress  tensor,  again  displaying  the  proper 
anisotropy  in  the  log  region  and  near  the  wall  and  the  lack  of  anisotropy  near  the  center- 
line.  Notice  also  that  the  normalized  shear  stress  r\2  and  the  wall-normal  Reynolds-stress 
component  do  behave  according  to  their  expected  asymptotic  behavior,  0(y)  and  0(y2) 
respectively. 

Figure  5.4(b)  shows  the  structure-dimensionality  components  for  this  case.  Of  note  is 
the  dli  component.  It  is  the  smaller  component  indicating  structures  preferentially  aligned 
with  the  ^-direction.  Furthermore  it  shows  a  minimum  near  the  wall,  where  near-wall 
streaks  aligned  with  the  flow  direction  have  come  to  be  expected. 

Figure  5.4(c)  shows  the  eddy-axis  tensor  components.  In  the  log  region,  in  locally  homo¬ 
geneous  shear  flow,  the  transverse  components  022  and  033  coincide  (as  in  plain  homogeneous 
shear  flow),  whereas  near  the  wall  0,22  is  forced  to  vanish  by  the  wall-blocking  applied. 

Finally  Figure  5.4(d)  displays  the  structure  parameters  used  in  the  model.  Note  that  for 
the  particular  case  of  no  frame  rotation  (3  is  modeled  as  unity  throughout  the  flow  region 
(in  this  shear  flow  situation),  and  that  0  grows  large  near  the  wall,  indicating  that  the 
near- wall  streaks  are  mostly  of  a  jetal  nature.  Thus  enforcing  <f>  =  1  at  the  wall  as  required 
by  the  asymptotic  behavior  of  the  Reynolds  stress  is  consistent  with  its  near-wall  behavior 
without  blocking. 


5.4  Channel  Flow  with  Spanwise  Frame  Rotation 

Here  fif  =  0  and  the  frame  rotation  is  aligned  with  the  mean  flow  rotation.  The 
ASBM  will  be  compared  against  the  DNS  of  Alvelius  (1999).  This  database  consists  of 
two  friction  Reynolds  numbers,  Rer  =  {180,360}.  For  Rer  =  360  there  are  results  for 
the  fixed  frame  and  for  a  rotation  number  Ro  =  0.22.  For  the  lower  Reynolds  num¬ 
ber,  Rer  =  180,  results  are  available  for  a  fixed  frame  and  for  several  rotation  numbers, 
Ro  =  {0.055,0.11,0.22,0.43,0.77}.  First  the  results  for  the  larger  Reynolds  number  are  in¬ 
vestigated.  Then  the  smaller  Rer ,  and  finally  skin  friction  coefficients  and  mass  flow  rates 
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are  compared.  Here  in  the  spanwise-rotating  frame  case  the  rotation  number  is  defined  as 


Ro  = 


Qf2h 
Ub  ’ 


(5.15) 


where  is  the  magnitude  of  the  frame  rotation  rate,  h  is  the  half  height  of  the  channel 
and  Ub  is  the  bulk  velocity  in  the  channel.  To  match  the  DNS  fig  <  0  was  used. 

Figure  5.5  shows  mean  velocity  profiles  normalized  by  the  bulk  velocity  at  Rer  =  360. 
The  fixed-frame  result  is  similar  to  that  reported  previously  with  the  database  of  Moser 
et  al  Adding  frame  rotation,  Ro  =  0.22,  the  velocity  profile  becomes  asymmetric  about 
the  centerline  of  the  channel.  In  spanwise  frame  rotation,  the  Coriolis  terms  drop  out  of 
the  mean  flow  equations  making  them  insensitive  to  direct  effects  of  frame  rotation.  The 
mean  flow  asymmetry  is  then  a  secondary  effect  due  to  the  effect  of  the  frame  rotation  on 
the  Reynolds  stresses  and  turbulence  structure.  The  ASBM  does  capture  the  asymmetry 
of  the  velocity  profile. 

In  Figure  5.6  the  ASBM  turbulence  intensities  are  compared  to  the  DNS  results.  Once 
again  the  ASBM  captures  the  anisotropy  of  the  Reynolds  stresses  and  its  dependence  on 
frame  rotation.  Notice  in  particular  the  fact  that  the  wall-normal  intensity,  v outgrows 
the  streamwise  intensity  v!  in  the  core  region  of  the  channel. 

Figure  5.7  shows  profiles  of  the  structure  dimensionality  tensor  for  the  rotating  case. 
Comparing  with  Figure  5.4(b)  it  is  clear  that  the  dimensionality  is  little  affected  by  the 
frame  rotation.  It  does  display  an  asymmetry,  but  this  results  directly  from  the  asymmetry 
in  the  mean  velocity  gradient.  There  are  no  dramatic  changes  as  in  the  Reynolds  stresses. 

The  lower  Reynolds  number  case,  Rer  =  180,  is  less  enticing  for  model  development 
because  it  may  still  be  contaminated  by  low  Reynolds  number  effects.  In  particular  the 
wall  region  is  felt  across  much  of  the  domain.  Nevertheless  it  is  rich  in  different  rotation 
numbers  available. 

Figure  5.8  displays  mean  velocity  profiles,  normalized  by  the  bulk  velocity,  for  several 
rotation  numbers,  at  ReT  =  180.  The  fixed  frame  result  is  again  similar  to  that  explored 
previously  with  the  database  of  Moser  et  al  With  increasing  frame  rotation  rate,  the 
velocity  profile  becomes  more  asymmetric  about  the  centerline  of  the  channel.  The  DNS 
displays  a  region  of  the  flow  with  a  nearly  linear  velocity  profile,  which  is  not  clearly  captured 
by  the  ASBM.  This  might  be  the  result  of  low  Reynolds  number  effects  creeping  into  the 
model  results. 

Figures  5.9(a-d)  present  the  turbulence  intensities  for  the  higher  rotation  numbers. 
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Along  with  the  anisotropy  of  the  stresses  it  is  interesting  to  observe  the  increasing  separation 
between  the  v '  and  the  uf  intensities  in  the  core  region  of  the  channel. 

The  structure  dimensionality  is  explored  for  the  same  four  rotation  numbers  in  Fig¬ 
ures  5.10(a-d).  It  is  still  remarkably  insensitive  to  the  rotation  rate,  except  for  the  larger 
rotation  rates  where  becomes  smaller,  indicating  elongated  structures  in  the  streamwise 
direction,  and  becomes  increasingly  larger,  indicating  also  a  growing  effect  of  flattening 
of  most  structures  in  the  mean  and  frame-rotation  directions. 

Figure  5.11  displays  wall  friction  velocities  for  all  the  Reynolds  numbers  and  rotation 
numbers  investigated  so  far,  normalized  by  the  wall  friction  velocity  for  the  fixed  frame 
case.  The  stable  suction  side  experiences  turbulence  decay,  as  evident  in  the  left  hand  side 
of  the  mean  velocity  and  turbulence  intensities  profiles  shown  earlier  (see  Figures  5.8  and 
5.9  for  example).  The  unstable  pressure  side  experiences  turbulence  growth  as  evident  in 
the  same  figures.  This  is  again  evident  here  in  the  form  of  lower  friction  velocities  in  the 
stable  side  and  larger  friction  velocities  in  the  unstable  side.  The  ASBM  predicts  a  more 
gradual  growth  of  the  skin  friction  asymmetry,  compared  to  the  DNS  which  seems  more 
abrupt. 

Figure  5.12  shows  the  bulk  Reynolds  number,  Ret  =  Ubh/v ,  where  E/&  is  the  bulk  ve¬ 
locity.  Evidently  there  is  a  mass  flow  increase  with  the  friction  Reynolds  number.  A  more 
subtle  effect  is  the  mass  flow  rate  growth  with  increasing  rotation  number.  For  the  higher 
Rer  ~  360,  the  ASBM  seems  to  capture  the  DNS  trend.  For  the  smaller  Reynolds  number, 
the  ASBM  seems  slightly  over-sensitive  to  the  rotation  number,  perhaps  due  to  excessive 
near  wall  influence  (low  Reynolds  number  effects). 

5.5  Channel  Flow  with  Streamwise  Frame  Rotation 

Here  =  0  and  the  frame  rotation  is  aligned  with  the  mean  flow  direction.  The  flow 
and  the  turbulence  respond  symmetrically  or  anti-symmetrically  to  the  direction  of  the 
rotation  applied  (anti-clockwise  or  clockwise),  depending  on  if  the  quantities  vanish  for  no 
frame  rotation  or  not.  The  ASBM  will  be  compared  against  the  DNS  of  Oberlack,  Cabot, 
and  Rogers  (1998).  This  database  consists  of  a  channel  flow  at  a  friction  Reynolds  number 
Rer  =  180  at  two  rotation  numbers,  Ro  =  {3.2,10}.  Here  in  the  case  of  a  streamwise- 
rotating  frame  the  rotation  number  is  defined  as 

&2h 
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where  Qf  is  the  magnitude  of  the  frame  rotation  rate,  h  is  the  half  height  of  the  channel 
and  uT  is  the  friction  velocity  in  the  channel.  In  the  DNS  fi{  >  0  was  explored. 

This  is  a  severe  test  for  the  ASBM  as  the  model  was  not  optimized  for  this  flow,  and 
furthermore,  the  rotation  is  not  in  the  plane  of  the  flow,  making  this  a  three-dimensional 
problem.  The  model  responds  to  rotation  in  a  way  determined  by  the  optimization  for  the 
case  of  spanwise  frame  rotation.  As  the  flow  physics  are  not  the  same  in  the  present  case, 
the  ASBM  will  not  perform  as  nicely,  as  already  predicted  when  analyzing  the  behavior  of 
the  algebraic  structure  equations  to  different  normalized  deformation  rates  of  streamwise 
frame  rotation  (see  Section  3.12.4).  Because  the  effects  of  frame  rotation  are  parameterized 
by  a  scalar,  it  cannot  distinguish  the  presence  of  a  different  mode  of  rotation. 

Figure  5.13(a)  shows  the  streamwise  mean  velocity  U  normalized  by  the  friction  velocity 
for  the  different  rotation  numbers.  The  DNS  shows  very  little  sensitivity  to  the  frame 
rotation.  In  fact  the  centerline  velocity  even  decreases  for  the  larger  rotation  rate.  In 
contrast  the  ASBM  displays  a  monotonic  growth  of  the  centerline  velocity  with  rotation 
number,  indicating  the  stabilizing  effect  the  rotation  has  on  the  model.  Smaller  turbulent 
levels  imply  less  friction  velocity  and  therefore  higher  centerline  velocities.  The  problem  is 
more  evident  for  the  higher  rotation  rate,  and  is  similar  to  the  one  found  by  Oberlack  et  al. 
exploring  differential  second  moment  models. 

Figure  5.13(b)  shows  the  spanwise  mean  velocity  W  resulting  from  the  frame  rotation 
rate.  The  presence  of  the  frame  rotation  excites  Reynolds  stress  components  that  normally 
do  not  exist.  These  components  in  turn  are  responsible  for  the  appearance  of  the  span- 
wise  velocity.  Again  as  in  the  differential  model  explored  in  Oberlack  et  al. ,  the  ASBM 
underpredicts  the  DNS  results. 

Figures  5.14(a)  and  5.14(b)  show  the  normal  Reynolds  stresses  for  the  two  different 
rotation  numbers.  The  stresses  from  the  DNS  are  hardly  affected  by  the  rotation,  while  a 
somewhat  increased  effect  can  be  seen  in  the  ASBM. 

Figures  5.14(c)  and  5.14(d)  show  the  shear  stresses  for  the  same  flows.  Here  the  DNS 
shear  stresses  R\$  and  R23  become  quite  energized  as  the  rotation  number  is  increased, 
while  the  ASBM  stresses  do  not  grow  at  the  same  rate.  Nevertheless  it  can  be  seen  that 
the  model  predicts  the  correct  sign  for  the  shear  stresses. 

Figures  5.15(a)  and  5.15(b)  show  the  structure-dimensionality  tensor  components  for 
the  rotation  rates  of  3.2  and  10,  respectively.  Again  the  dimensionality  tensor  is  relatively 
insensitive  to  the  frame  rotation.  Nonetheless  a  reduction  of  is  observed,  indicating 
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an  elongation  of  the  structures  in  the  streamwise  direction.  The  off-plane  dimensionality 
components  d J3  and  have  been  activated,  but  they  are  just  barely  evident. 


5.6  Summary  and  Model  Assessment 

This  chapter  explored  the  complete  algebraic  structure-based  model,  consisting  of  the 
algebraic  equations  for  the  turbulence  structure  laid  down  in  Chapter  3  were  coupled  to  the 
scale  evolution  equations  from  Chapter  4.  The  chapter  concentrated  on  the  exploration  of 
wall-bounded  channel  flows,  under  different  imposed  modes  of  frame  rotation. 

The  homogeneous  structure  based  models  have  proved  capable  of  giving  realistic  repre¬ 
sentations  for  the  Reynolds  stresses  and  structure  dimensionality  tensors  for  different  mean 
deformations  with  or  without  frame  rotation  (see  Chapter  3,  and/or  previous  algebraic 
structure  based  models,  Reynolds  et  al  2000  and  Haire  2003). 

This  was  the  first  time  the  full  model  was  implemented  in  a  wall-bounded  flow,  integrat¬ 
ing  all  model  equations  through  the  viscous  sublayer  (say  <  30).  A  set  of  fixed  frame 
channel  flows  was  explored  to  investigate  the  model  performance  when  dealing  with  near 
wall  regions.  The  results  are  encouraging  as  good  log  laws  were  achieved,  as  well  as  proper 
anisotropies  for  the  Reynolds  stress  and  dimensionality  tensors  in  the  log  region.  However 
near  the  centerline,  when  the  mean  deformation  vanishes,  the  model  fails  to  predict  any 
anisotropy,  contrary  to  DNS  results.  This  is  a  feature  common  to  all  algebraic  models, 
reflecting  the  local  nature  of  their  formulation. 

Model  predictions  for  the  spanwise-rotating  frame  again  did  not  show  any  spurious  effect 
of  the  integration  through  either  viscous  sublayer:  the  stable  side  can  be  quite  challenging, 
with  its  small  levels  of  turbulence  activity.  The  model  did  display  asymmetric  mean  velocity 
profiles  and  proper  anisotropy  levels,  including  the  cross  over  of  the  normal  stresses.  The 
qualitative  comparison  with  DNS  reflects  the  adequacy  of  the  choices  for  the  structure 
parameters  made  in  Chapter  3. 

The  model’s  response  to  increased  frame  rotation  rate  in  terms  of  wall  skin  friction  is 
a  little  too  gradual  at  first,  but  then  saturates  as  in  the  DNS.  The  mass-flow  rate  grows  a 
little  faster  than  the  DNS.  It  is  not  clear  however  what  effect  the  low  Reynolds  number  of 
the  simulations  might  have  on  these  results. 

Model  predictions  for  the  streamwise-rotating  frame  again  show  adequate  behavior  as 
the  viscous  sublayer  is  bridged.  This  time  however  the  model  did  display  an  oversensitivity 
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to  the  rate  of  frame  rotation.  A  stabilizing  effect  of  rotation  on  the  turbulence  is  observed 
in  the  model  results,  opposite  to  the  rather  insensitive  character  of  the  DNS  results.  This 
fact  is  blamed  in  the  fact  that  the  model  (i)  was  optimized  for  the  case  of  spanwise  frame 
rotation,  which  required  a  rather  large  parameter  sensitivity  to  frame  rotation,  and  (ii)  the 
scalar  parameterization  chosen  cannot  transmit  information  about  the  direction  of  frame 
rotation  back  to  the  model  for  obvious  reasons.  Also  surprising  is  the  fact  that  the  secondary 
flow  generated  is  much  weaker  than  predicted  by  DNS,  perhaps  a  casualty  of  the  excessive 
stabilization  driven  by  the  frame  rotation. 
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Figure  5.1:  Equilibrium  ratio  of  production  to  dissipation  of  turbulent  kinetic  energy  for  homoge¬ 
neous  shear  in  a  rotating  frame 


Figure  5.2:  Mean  velocity  for  the  fixed-frame  channel  at  Rer  =  {180,395,590,3000}.  Solid  lines: 
DNS.  Dotted  line:  logy+/0.41 -1- 5.  Dashed  line:  y+  and  log2/+/0.38  +  4.1.  Dot-dashed  line:  ASBM. 


Figure  5.3:  Turbulence  intensities  and  shear  stress  for  fixed-frame  channel  at  Rer  =  590.  Symbols: 
DNS  (Moser  et  al.  1999).  Solid  line:  ASBM. 
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Figure  5.4:  Turbulence  structure  and  model  pai 
Symbols:  DNS  (Kim  1992).  Solid  lines:  ASBM. 
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Figure  5.5:  Mean-velocity  normalized  by  bulk-velocity  for  the  spanwise-rotating  channel  at  Rer  — 
360  and  Ro  =  {0,0.22}.  Symbols:  DNS  (Alvelius  1999).  Solid  lines:  ASBM. 
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Figure  5.6:  Turbulence  intensities  and  shear  stress  for  the  spanwise-rotating  channel  at  Rer  =  360, 
and  Ro  =  0.22.  Symbols:  DNS  (Alvelius  1999).  Solid  lines:  ASBM. 
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Figure  5.7:  Turbulence  structure-dimensionality  for  the  spanwise-rotating  channel  at  Rer  =  360  and 
Ro  =  0.22.  Solid  lines:  ASBM. 
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(a)  Ro  =  0.11 


(b)  Ro  =  0.22 


(c)  Ro  =  0.43  (d)  Ro  =  0.77 


Figure  5.9:  Turbulence  intensities  and  shear  stress  for  spanwise-rotating  channel  at  ReT  =  180  and 
Ro=  {0.11,0.22,0.43,0.77}.  Symbols:  DNS  (Alvelius  1999).  Solid  lines:  ASBM. 
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(a)  Ho  =  0.11 


(b)  Ho  =  0.22 


/  / 
(c)  Ho  =  0.43  (d)  Ho  =  0.77 


Figure  5.10:  Turbulence  structure-dimensionality  for  the  spanwise-rotating  channel  at  ReT  =  180 
and  Ro  =  {0.11,0.22,0.43,0.77}.  Solid  lines:  ASBM. 
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Ro 


Figure  5.11:  Friction  velocity  vs.  rotation  number.  *,  +:  DNS  (Alvelius  1999).  □,  ■:  ASBM. 
Spanwise-rotating  channel  at  ReT  =  {180,360}  respectively. 
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Figure  5.12:  Bulk-velocity  Reynolds-number  vs.  rotation  number.  *:  DNS  (Alvelius  1999). 
□:  ASBM.  Spanwise-rotating  channel  at  ReT  =  {180,360}  respectively. 
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Figure  5.13:  Mean  streamwise- velocity,  f/+,  and  spanwise- velocity,  W+,  for  streamwise- rotating 
channel  at  Rer  —  180  and  Ro  —  {3.2, 10}.  Symbols:  DNS  (Oberlack  et  al.  1998).  Lines:  ASBM. 
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Figure  5.14:  Reynolds  stresses  for  streamwise-rotating  channel  at  Rer  =  180  and  Ro  =  {3.2, 10}. 
Symbols:  DNS  (Oberlack  et  al  1998).  Lines:  ASBM. 


Chapter  6 


Conclusions  and  Recommendations 


6.1  Summary  and  Conclusions 

The  starting  point  of  this  work  is  the  appreciation  that  there  is  more  than  one  indepen¬ 
dent  tensor  variable  describing  a  turbulence  state.  For  instance,  the  Reynolds  stress  tensor 
and  the  structure  dimensionality  tensor  carry  independent  information.  With  mean  and/or 
frame  rotation  both  are  important.  Structure  information  is  brought  into  the  proposed 
turbulent  model  through  the  use  of  hypothetical  turbulent  eddies,  each  assigned  a  set  of 
properties.  Averaging  over  an  ensemble  of  eddies  produces  eddy  statistics,  and  constitutive 
relations  between  these  statistics  and  the  Reynolds  stress  and  dimensionality  tensors. 

The  proposed  model  uses  an  algebraic  approximation  for  the  eddy  statistics.  This 
implicitly  assumes  the  existence  of  an  analysis  frame  where  the  turbulence  is  locally  in 
(or  near)  equilibrium  with  the  mean  flow.  The  farther  away  from  equilibrium,  the  worse 
the  algebraic  approximation  will  be.  The  algebraic  model  produces  physically  consistent 
turbulence  states  for  different  combinations  of  normalized  mean  strain  and  mean  or  frame 
rotation  rates,  and  displays  the  property  of  material  indifference  of  2D  turbulence  to  rotation 
aligned  with  its  direction  of  independence.  A  weakness  appears  for  the  case  of  shear  in 
presence  of  streamwise  frame  rotation.  Here  the  principal  shear  stress  is  too  sensitive  to 
the  frame  rotation,  while  the  secondary  stresses  are  underpredicted. 

The  algebraic  equations  for  the  turbulence  structure  are  complemented  by  a  set  of 
scale  equations  for  the  turbulent  kinetic  energy  and  the  large-scale  enstrophy,  originally 
developed  for  homogeneous  flows  by  Reynolds,  Langer,  and  Kassinos  (2002).  Exact  terms 
in  the  evolution  equations  represent  large  scale  processes,  leaving  the  modeled  equations 
with  a  more  solid  basis  than  traditional  choices,  especially  for  the  second  scale. 

Both  the  algebraic  equations  and  the  scale  equations  are  sensitized  to  the  presence  of 
a  wall  in  a  way  that  preserves  proper  asymptotic  behavior  for  different  quantities  in  the 
vicinity  of  a  wall.  The  essence  of  the  idea  is  to  get  the  relative  magnitude  of  different 
quantities  properly.  For  example,  in  the  linear  sublayer  (say  y+  <  5)  molecular  transport 
should  dominate  over  turbulent  transport. 
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The  complete  model  —  evolution  equations  for  the  turbulent  scales,  and  algebraic  equa¬ 
tions  for  the  energy-normalized  turbulence  structure  —  was  investigated  and  optimized  for 
a  set  of  channel  flows  in  fixed  and  rotating  frames  of  reference.  This  is  the  first  time  a  com¬ 
plete  algebraic  structure  based  model  has  been  applied  to  the  solution  of  a  wall  bounded 
flow.  For  fully  developed  channel  flow,  using  a  single  set  of  model  coefficients,  it  was  found 
that: 

•  the  near- wall  modifications  successfully  allowed  the  integration  of  the  model  equations 
through  the  viscous  sublayer  (say  y+  <  30)  without  indication  of  any  spurious  effects; 

•  in  the  fixed-frame  case,  the  mean  velocity  and  the  Reynolds  stresses  compare  favorably 
with  DNS,  near  the  wall  and  in  the  logarithmic  region; 

•  the  model  responded  favorably  to  spanwise  frame  rotation,  displaying  an  increased 
asymmetry  of  the  velocity  profile  with  rotation  number,  similar  to  the  DNS  results. 
The  Reynolds  stresses  are  roughly  in  agreement  with  the  DNS,  with  particular  mention 
to  the  crossover  of  the  normal  stresses  in  the  unstable  side  of  the  channel; 

•  in  a  streamwise-rotating  frame,  the  model  displayed  excessive  turbulent  stabilization 
for  large  rotation  numbers  while  underpredicting  the  secondary  motion  resulting  from 
the  frame  rotation  (reminiscent  of  second  moment  model  results  reported  along  with 
the  DNS  results), 

•  in  all  flows,  when  the  mean  deformation  vanishes,  the  model  predicts  an  isotropic 
state,  in  contrast  to  DNS  results,  due  to  the  inability  of  the  algebraic  model  to  deal 
with  non-equilibrium  situations. 

6.2  Recommendations  for  Future  Work 

This  work  primarily  addressed  model  development  for  diverse  deformation  rates  with 
and  without  frame  rotation,  and  model  development  for  compliance  with  near-wall  regions. 
There  remain  areas  that  need  further  exploratory  work. 

•  The  algebraic  model  for  the  eddy-axis  tensor  needs  to  be  revisited  for  compliance  with 
required  behavior  near  the  edge  of  a  shear  layer  (turbulent /non-turbulent  interface), 
as  suggested  in  Equation  (4.64)  for  example. 
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•  The  model  needs  to  be  assessed  in  increasingly  complex  flows,  culminating  with  the 
evaluation  of  its  performance  in  rotating  machinery  configurations. 

•  Explore  some  analytical  approximation  to  incorporate  non-local  effects  when  the  mean 
deformation  vanishes. 

•  Revisit  the  algebraic  equations  to  address  the  excessive  turbulence  stabilization  that 
develops  for  shear  in  streamwise-rotating  frame. 

When  the  structure  approach  is  deemed  robust  and  tractable  enough  for  complex  engi¬ 
neering  flows,  a  structure-based  extension  for  the  energy  equation  should  be  pursued.  This 
might  for  example  involve  assigning  an  additional  temperature  property  to  the  individual 
eddies,  in  the  same  way  other  eddy  properties  were  assigned  in  Section  3.3,  and  then  ana¬ 
lyzing  how  averaged  eddy  properties  combine  to  produce  structure-based  statistics,  useful 
for  the  averaged  mean  flow  and  energy  equations. 


Appendix  A 


Flattened  Eddy- Axis  Model 

This  appendix  details  the  construction  of  the  structure  model,  following  Reynolds  (personal 
communication). 

In  the  eddy-axis  model,  a  turbulence  field  is  described  in  terms  of  an  ensemble  of  hy¬ 
pothetical  eddies.  Each  eddy  is  described  in  terms  of  an  eddy-axis  unit  vector,  #,  and  an 
eddy-flattening  unit  vector,  6,  that  gives  a  non-axisymmetric  form  to  the  eddy.  Given  an 
eddy-axis  vector  and  a  flattening  direction,  the  turbulence  structure  is  expressed  in  terms 
of  moments  conditioned  on  both  vectors,  and  this  is  used  as  basis  for  model  statistics.  The 
first  step  is  developing  an  expression  for  the  conditioned  dimensionality  tensor.  This  is  then 
later  used  in  the  development  of  the  conditioned  Reynolds-stress  tensor. 

A.l  Conditioned  Dimensionality  Tensor 

The  conditioned  dimensionality  tensor  represents  an  average  over  all  eddies  that  share 
the  same  direction,  with  same  flattening  orientation.  It  is  first  written  as  the  most  general 
function  of  the  two  vectors,  a  and  b.  Appropriate  constraints  will  then  determine  the  free 
coefficients. 

The  most  general  form,  given  two  vectors,  is 

(V^UiTijlS,  b) /  V2  =  C\ Sij  -f-  C2Q'i&j  ■+■  Czbibj  4-  C/^{pL%bj  +  djbi)  (3*3) 

+  C^{€ipq€ipbqbj  +  CjpqCLpbqbi)  -|-  C^i^ipqdpbqdj  4*  tjpq&pbq&i) , 

where  V2  —  (F2|a,6). 

The  orthogonality  between  the  eddy-axis  tensor  and  the  normal  vector  is  expressed  by 
the  constraint 

djUj  =  0.  (A.8) 

Multiplying  Equation  (3.3)  by  dj  gives  relations  between  constants  (recall  djdj  =  1) 
dj(V2nirij\d, b)/  V2  =  0  =  C\di  +  C^di  +  Cz(dnbn)bi  +  C4[(dnbn)di  +  6J  (A.9) 

“h  G5  (djibji)€ipqdpbq  ~h  C^ipqdpbq. 
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Equating  coefficients  of  like  terms  gives 

Cl  -f-  C*2  +  ( anbn)C 4  =  0,  (dn^n)Cs  +  C4  =  0,  (fln^n)C^5  +  Cej  =  0, 


or 


C6  = 

{Q'nbn^Cs,  C4  —  {p'n^n)C^1  C2  —  C\  (Un^n)  C$, 

(A.  10) 

The  next  constraint  is  the  unit  trace  (normalization)  requirement, 

<1 V2nm\a,  b)/  V2  =  1  =  3Ci  +  C2  +  C3  +  2 C4{anbn). 

(All) 

Substituting  Equations  (A.  10)  the  normalization  condition  becomes 

20!  +  [1  -  (< %„bn)2}C3  =  1. 

(A.  12) 

Defining 

II! 

O 

III 

(A.  13) 

the  other  constants  can  be  represented  by 

Cl  =  H1  —  x[l  -  (anbn)2}}, 

(A.  14) 

C2  =  — 1{1  -  X[1  +  (Mn)2]}, 

(A.  15) 

C4  —  ~-{Q'nbn)X'> 

(A.16) 

Cq  ~  —  (Un^n)C* 

(A.17) 

Hence 

{V2ni7ij\d 

i,b}/  V 2  =  ±{1  -  x[l  -  [anbn)2}}8ij  -  ±{1  -  x[l  +  (anbn)2}}aiaj 

(A.18) 

+  xkbj  -  ( anbn)x(aibj  +  ajbi) 

Q{tipqQ>pbqbj  H“  CjpqQ'pbqbi)  ( ^rJ^njCi^ipq^p^q^j  4”  ^-jpq^p^q^i) • 


The  most  basic  eddy-axis  model  is  an  axisymmetric  model  which  does  not  include  any 
flattening  effects,  i.e.  x  =  C  =  0.  It  seems  only  natural  to  start  analyzing  flattening  by 
adding  one  perturbation  at  a  time.  That  is  achieved  by  allowing  a  nonzero  x  and  taking 

C  =  o, 


{V2ninj\a,b)/  V2  =  ±{1  -  *[1  -  (an&„)2]}<%  -  §{1  -  *[1  +  (a^)2]}^- 
+  Xpibj  —  (anbn)x(aibj  +  o,jbi). 


(3.4) 
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A.l.l  Special  Cases 

If  X  =  0  the  axisymmetric  limit  is  obtained, 

(V2rii7ij\a,  b)/  V2  =  —  didj).  (A.19) 

If  dnbn  =  0,  the  eddy-axes  are  all  perpendicular  to  the  vector  that  organizes  the  flatten¬ 
ing,  and 

{V2mnj\a,  b)f  V2  —  |(1  —  x)(Sij  -  a^j)  +  xhbj. 

Furthermore,  if  \  =  1,  the  eddies  are  fully  flattened  in  the  direction  of  the  b  vector, 

{V2nirij\d,b) /  V2  =  bfij. 

A.  2  Conditioned  Reynolds- Stress  Tensor 

Following  Kassinos  and  Reynolds  (1994),  hereafter  referred  to  as  TF-61,  the  velocity 
vector  is  written  in  an  eddy  coordinate  system,  with  a  component  normal  to  the  eddy  and 
a  component  aligned  with  the  eddy-axis  direction, 

—  V{vmdrn)di  4"  V (trst'Vr'R'sdtjtipqdpTlq'  (3*5) 

Taking  the  second  moment  of  the  velocity  gives  an  expression  for  the  conditioned  stresses 
in  terms  of  the  eddy-axis  vector  and  the  conditioned  dimensionality  tensor, 

=  {ViVj\a,  b)  -  ((vmam)2V2aiaj |o,  b)  (3.6) 

4*  (  tipqdpd  jTiq  4-  ejpqdpdinq)  |  a,  b) 

4"  ((€rstVrrl8Q>t)^V‘*£ipq€jklQ'p'n'qG'knl\Q‘)  &)• 

To  simplify  Equation  (3.6),  the  identity 

^rst^opq  =  ^ro^sp^tq  +  Srp^sq^to  4"  firq^so^tp  8ro5SqStp  —  Srp^so^tq  “  ^rq^sp^to  (2.11) 
is  used.  Then  the  coefficient  of  the  third  term  in  Equation  (3.6)  simplifies  to 


(erS£'yr7i5atf)2  —  ( €rstVrnsdt)(€0pqV0npdq ), 
—  1  (^m®m)  • 


(A.20) 
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Next  (j)  and  %  are  defined  as 

0  = 

'ygTlg  =  (€r3tVr(lt){vTnQ/Tn)Tl3 
Finally,  part  of  the  third  term  in  Equation  (3.6)  can  be  simplified  as 

tipqtjklQ'pTlqQ'k'H'l  Cli&j  7l%7lj.  (A.21) 

Substituting  Equations  (A.20),  (3.7),  (3.8),  and  (A.21)  into  Equation  (3.6)  gives 

Rf/  =  (j>V2  aia j\a,b)  (A.22) 

“I"  HsV  (EipqQpdjTlqTlg  H-  €jpqdpd%TtqTlg)\d^b') 

+  {(1  -  <i>)V2(6ij  -  djdj  -  UiUj) |a,  b). 

Substituting  the  conditioned  dimensionality  tensor,  Equation  (3.4),  into  Equation  (A.22) 
produces 

Rf/  =  (V24>diaj)  +  (V2(l  -  m(6ij  -  aiaj)  +  x{J[l  -  (anbnfjSij  (3.9) 

-  i[l  +  (dnbn)2]didj  -  bibj  +  {dnbn)(dibj  +  djbi)}]) 

+  {V2'y3(6ipqdpdj  +  6jpqdpdi)\^  { 1  —  x[l  (dnbn)  ]}^s  +  Xbqbg  x(®ra^n)(®s^<j)])- 

The  final  step  in  the  modeling  process  requires  taking  the  ensemble  average  over  all 
eddies  (eddy-axis  orientations  and  flattening  directions,  a  and  b).  TF-61  introduces  the 
definitions 

q2rij  =  ( ViVj ),  q2dij  =  ( V2mnj ),  q2dij  =  {V2didj),  q\j  =  ( V2bibj ),  (3.10) 

q2<j>  =  ( V2<t> ),  q2jk  =  (V2%),  q2x  =  (’ V2X ),  (3-H) 

respectively  the  energy- weighted  normalized  Reynolds  stress  tensor,  eddy-axis  tensor,  eddy¬ 
flattening  tensor,  and  the  energy- weighted  structure  jetal  parameter,  structure  helix  vector, 
and  structure  flattening  parameter. 

Higher  moments  also  appear  during  the  averaging  of  Equation  (3.9).  Among  others, 
(1 V2<j)aiaj )  and  (V27 8eipqapajx{anbn)(asbq)) .  These  are  modeled  as  if  the  scalars  were  inde¬ 
pendent  from  the  tensorial  forms,  and  with  special  care  to  keep  some  of  the  traces  correct. 


(3.7) 

(3.8) 
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In  the  two  examples  mentioned  the  results  are 

(V24>ai(ij)  =  q2(j)aij , 

(V  ^s^ipqdpdjX{dnbn){dsbq))  =  (V  ^s^ipqdpdjX{dsO'nbnbq)) 5 

=  #27s  £ipq  apj  X  {asjfinq)' 

The  approximation  is  exact  for  some  limiting  cases.  In  irrotational  flow  in  a  fixed  frame  of 
reference,  (f>  =  7*  =  0,  and  both  sides  are  zero.  The  case  of  no  mean  flow  deformation  is  a 
subset  of  the  one  just  mentioned.  On  the  other  side  of  the  spectrum,  in  rapidly  deformed 
flows  in  which  the  eddy-axis  tensor  approaches  the  1C  state  the  approximations  again  apply. 

With  this  approximation,  and  definitions  given  by  Equations  (3.10)  and  (3.11),  the 
ensemble  averaged  normalized  structure  tensors  become 

dij  “  X(1  —  dnmbmn)\5ij  (3.12) 

—  —  X(1  4  anmbmn)\^ij  +  X^ij  ~  X(ainbnj  4-  djnbn%\ 

rij  =(1  —  0)  2  i^ij  ~  aij )  4  faij  (3.13) 

4(1““  <j>)X  [|(1  ””  anm^mn)^ij  ~  2  (1  4  dnrnbmn)aij  —  ft ij  4  din^nj  4  djnbnj\ 

4  (—7 /QT)(€iprdpj  4  £jprdpi){\[  1  ~  X(1  “  anm^mn))^kr  4  X^kr  ~~  Xakn^nr}- 

fij  is  not  an  independent  quantity.  It  can  be  obtained  through  the  relation 

=  Sij . 


rij  4  dij  4  fij 


(2.19) 
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Algebraic  Model  for  the  Eddy- Axis  Tensor 


This  appendix  details  basic  thoughts  involved  in  the  construction  of  the  model  for  the  eddy- 
axis  tensor,  as  it  appears  in  Reynolds,  Kassinos,  Langer,  and  Haire  (2000).  Initial  attempts 
at  writing  down  a  model  for  as  a  single  expression  in  terms  of  the  mean  strain  and 
mean  rotation  always  showed  the  undesirable  feature  of  producing  unrealizable  eddy-axis 
tensor  states  for  some  combinations  of  mean  strain  and  mean  rotation.  To  overcome  this 
issue,  the  computation  of  the  eddy-axis  tensor  was  separated  in  two  parts.  First  a  model 
that  takes  into  consideration  only  the  mean  strain  of  the  flow  is  proposed.  This  form  is 
found  to  produce  realizable  states.  In  a  second  step  the  eddy-axis  tensor  is  sensitized  to 
rotation  through  a  rotation  operation,  function  of  the  mean  rotation,  which  preserves  the 
realizability  characteristic  of  the  initially  strained  eddy-axis  tensor. 


B.l  Eddy- Axis  Tensor  Under  Pure  Strain 


The  first  part  of  the  model  deals  only  with  pure  strain,  using  the  simple  eddy-axis  model 

1  r  ,  (S*kakj  +  Sjkaki  —  §S 'mnanmSij)T 

Oij  —  «  ______  . 

3  ao  +  ai  y/r2S^pS*kqapq 

In  the  principal  coordinates  of  the  anisotropic  strain-rate  tensor  Stj  =  Sij  -  SkkSij/Z,  where 
Sij  —  $(Gij  +  Gji) ,  take 

t  r  o  o  ^ 

0  -i(T  +  o)  0 
\0  0  -i(T-o)  ) 

For  T  >  0  and  a  <  2,  St  1  dominates  and  in  the  RDT  limit  when  Tr  »  1, 


c*  _ _ 

Di j  ”  2 


aij  — 


/  1  o  o  ^ 
0  0  0 
0  0  0 


Snrrfimn  —  F- 
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B.2  Effect  of  Mean  Rotation 

Suppose  there  is  a  combination  of  strain  and  rotation.  The  eddy-axis  tensor  is  first 
calculated  just  considering  strain,  denoted  by  afj. 

(3.15) 

In  order  to  generate  a  realizable  a^-  including  effects  of  rotation,  a  rotation  operation  is 
performed  on  a\p 


s  lr  .  +  S*kal 

aij  =  ody  + - r= 

3  ao  +  2./r^Sj 


3  mn  ®nm  ^ ij  ) 


'*  S*  ns 
kp°kqapq 


&ij  —  HikHjl&ki- 


(3.16) 
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where  the  rotation  tensor  Hij  must  satisfy  the  orthonormal  conditions 

HikHjk  —  Sjj,  HkiHkj  =  Sij • 


is  modeled  as  a  function  of  the  mean  rotation-rate  tensor  fijj, 


Hij  —  C\  S{j  -}-  h\  C^Qij  “I”  C3 ^ik^kj ' 


The  coefficients  are  determined  by  the  orthonormal  condition, 

HikHjk  =  =  CfSij  +  2CiC3n?.  -  Cfflj  +  C§iknf%j, 

where  Qfj  =  flik^kj-  Now  a  Cayley-Hamilton  result  is  needed: 


^12  =  ^12^21^12  +  ^12^23^32  +  ^13^31^12  =  “^12(^12  +  ^13  +  ^23)  =  —  ~^2^12? 


where  Q2  =  Slnk£lnk-  Also, 


fifl  =  ^12^23^31  +  ^13^32^21  =  0. 


Hence  the  Cayley-Hamilton  result  is  fi?-  =  —  ^n2^-,  and  the  orthonormality  condition 
gives 


HikHjk  =  Sij  =  C\Sij  +  2C1C3n?-  -  Cl&ij  -  \C\tfVL%. 
Matching  the  tensors, 

Ci  =  i,  2CiC3  -cl-  ic|n2  =  o, 


Letting 


Cl  =  Cfft2/ 2. 


c2  =  ±-n  x/2/3  - 1/?. 


then 
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The  parameter  (3  will  be  set  by  reference  to  RDT  for  combined  plane  strain  and  rotation, 
where 


(  0 

e 

0  ^ 

0 

UJ 

0  ^ 

e 

0 

0 

p 

II 

—UJ 

0 

0 

1° 

0 

0  J 

0 

0 

0  J 

The  equations  for  a  material  line  element  L  in  the  mean  flow  are 


L\  —  G12L2  —  (e  +  ct>)Z/2,  T2  —  G21L1  —  (e  —  ui)L\. 


These  produce 


(e  +  uj)  - 


(e-a;). 


The  fixed  point  solution  gives 


an 


L 2  e  H-  uj  L2  e  —  ctJ  L\L2 

L*  =  ~2e~'  °22  =  I5  =  ~~2e~’  °12  = 


\/e2  —  a;2 
2e 


for  (e)  —  a;  >  0. 


For  this  case 


n2  =  2w2,  C3  =  ^,  C2  =  ±lv//?-W, 


so  the  rotation  matrix  is 


Hn  =  H22  =  1  -  /3/2,  #12  =  -H2l  =  ±^p-p/4. 


Then 


an  =  HlkHuaskl  =  i  ±  (1  -  0/2W0-  Wfi, 


022  =  J=F(1-/3/2)v'/3-/32/4, 


012  =  4  -  P  +  /?2/4. 

Note  that  on  +  022  =  1  for  any  (3 .  Let  r  =  uj2/e2.  Then  from  the  RDT  limit  solution 

012  =  |\/l -r  =\-P  +  P2/  4. 
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Solving  for  6. 

/3  =  2±2y/  \  +  i-v/T^r. 

The  minus  sign  is  taken  in  order  that  /3  =  0  when  r  =  0  (no  rotation), 

p  =  2-2^JT+Jv^^. 

an  -  G22  gives 

an  —  &22  =  ±\/r  =  ±2(1  —  f3/2)y/ (3  -  /32/4. 

Using  the  (3  solution, 

l-/?/2  =  ^[1  +  a/T^], 

/3-/32/4=i[l-vT^]. 

Hence 

an  ~  ^22  =  =  ±2^/|  [1  +  Vl  -V]  [1  -  VT^i)  =  ±Vr, 

so  the  analysis  checks.  Hence  for  r  <  1 

P  =  2  —  2^\\l  +  Vn^\,  r  <  1. 

The  sign  ambiguity  in  C2  is  resolved  by  looking  at  an 

an  =  ${l  +  r)  =  i±  (1  -  0/2)  V/?  -  W*  =  k  ±  Jr, 
so  the  +  sign  is  required.  The  two  coefficients  are  therefore 

cs=|. 

Now  a  generalization  is  made  replacing  r  by 

apqQqrSrp 

^kn^nmamk 

In  the  RDT  limit  used  above,  this  definition  gives  the  same  value  for  r  used  above,  r  =  u2/ e2. 
This  definition  of  r  allows  the  model  to  be  extended  to  elliptic  streamline  flows,  where 
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\co\  >  e ,  and  to  rotations  in  other  planes.  It  also  produces  material  indifference  when  the 
rotation  is  about  an  axis  of  independence.  Apparently  r  >  01. 

Consider  what  happens  when  J1  »  e2  (elliptic  mean  streamlines).  If  f3  =  2  is  chosen 
for  large  r,  then  Hn  =  H2 2  =  0,  #21  =  1,  and  Hu  =  -1,  which  produces  an  —a22~\  and 
a12  =  This  says  that  the  dominant  (most  energetic  eddies)  are  at  right  angles  to  the 
axis  of  principal  strain  rate,  which  at  first  glance  may  seem  unreasonable.  However,  the 
rate  of  energy  input  to  the  eddy  varies  as  cos  (a>£),  where  t  —  0  corresponds  to  the  time  at 
which  the  eddy  axis  is  aligned  with  the  axis  of  principal  strain  rate.  Hence 

dk 

—  =  AcosUvt)  k  —  Bsm(cjt) 
dt 


so  it  is  indeed  reasonable  that  in  these  strongly  elliptic  flows  the  most  energetic  eddies  lag 
the  angle  of  maximum  straining  of  the  eddies  by  7r/2. 

Guided  by  the  discussion  above  the  following  generalization  is  proposed 


/?  = 


2  -  2^/|(l  +  vT^T)  if  r  <  1, 


Using  a  =  y/2j3  -  /32/2, 


a=  < 


2  —  2^j i(l  —  i/l  —  1/r)  if  r  >  1. 


\/l  —  y/1^7  if  r  <  1, 
yj  1  +  y/l  —  l Jr  if  r  >  1. 


The  model  described  above  is  used  to  generate  the  eddy-axis  tensor,  a^-,  in  the  frame  for 
which  the  mean  deformation-rate  is  constant  The  argument  for  this  is  that  the  evolution 
of  the  eddy-axis  tensor  is  kinematic,  and  since  the  algebraic  model  is  an  equilibrium  model 
it  should  be  applied  in  the  frame  for  which  the  mean  deformation  rate  is  essentially  steady. 
Hence  for  shear  in  a  rotating  frame  (turbomachinery)  we  apply  the  model  above  in  the 
frame  of  the  rotor  or  stator.  However,  in  the  Reynolds  stress  model  the  parameters  are 
based  on  the  total  rotation  rate  (frame  rotation  rate  plus  mean  rotation  rate  relative  to  the 
frame)  because  the  dynamics  of  the  velocity  in  a  given  eddy  is  determined  by  the  absolute 
rotation  rate. 

XWCR  has  found  r  >  0  for  all  cases.  However,  during  the  iteration  for  the  solution  r  <  0  can  be 
encountered.  The  program  uses  a  smooth  continuation  of  the  functions  to  negative  r  in  order  to  achieve  a 
convergent  iteration,  and  returns  an  error  flag  if  the  final  r  is  negative 


References 


Alvelius,  K.  (1999).  Studies  of  Turbulence  and  its  Modelling  through  Large  Eddy -  and  Di¬ 
rect  Numerical  Simulation .  Ph.  D.  thesis,  Department  of  Mechanics,  KTH,  Stockholm, 
Sweden. 

Bernard,  P.  S.  (1990).  Turbulent  vorticity  transport  in  three  dimensions.  Theor.  and 
Comp .  Fluid  Dyn.  2 ,  165-183. 

Blaisdell,  G.  A.  and  K.  Shariff  (1996).  Simulation  and  modeling  of  the  elliptic  streamline 
flow.  In  Studying  Turbulence  Using  Numerical  Simulation  Databases  VI.  Proceedings 
of  the  1996  Summer  Program ,  pp.  433-446.  NASA  Ames/Stanford  Center  for  Turbu¬ 
lence  Research. 

Cazalbou,  J.  B.  and  P.  Bradshaw  (1993).  Turbulent  transport  in  wall-bounded  flows. 
Evaluation  of  model  coefficients  using  direct  numerical  simulation.  Phys.  Fluids  5(12), 
3233-3239. 

Cazalbou,  J.  B.  and  P.  Chassaing  (2002).  The  structure  of  the  solution  obtained  with 
Reynolds-stress-transport  models  at  the  free-stream  edges  of  turbulent  flows.  Phys. 
Fluids  14(2),  597-611. 

Cazalbou,  J.  B.,  P.  R.  Spalart,  and  P.  Bradshaw  (1994).  On  the  behavior  of  two-equation 
models  at  the  edge  of  a  turbulent  region.  Phys .  Fluids  5(5),  1797-1804. 

Coles,  D.  E.  and  E.  A.  Hirst  (1969).  Memorandum  on  data  selection.  In  D.  E.  Coles 
and  E.  A.  Hirst  (Eds.),  Proceedings  of  the  1968  AFOSR-IFP-Stanford  Conference  on 
Computation  of  Turbulent  Boundary  Layers ,  Volume  II,  Stanford,  CA,  pp.  47-54. 

Durbin,  P.  A.  and  B.  A.  Pettersson  Reif  (2002).  The  elliptic  relaxation  method.  In 
B.  Launder  and  N.  Sandham  (Eds.),  Closure  Strategies  for  Turbulent  and  Transi¬ 
tional  Flows ,  pp.  127-152.  Cambridge  University  Press. 

Greenspan,  H.  P.  (1968).  The  Theory  of  Rotating  Fluids.  Cambridge  University  Press. 

Haire,  S.  L.  (2003).  Toward  An  Affordable  Two-Equation,  Structure- Based  Turbulence 
Model  Ph.  D.  thesis,  Mechanical  Engineering  Department,  Stanford  University. 

Jeffreys,  H.  (1931).  Cartesian  Tensors.  Cambridge  University  Press. 


115 


116 


Kassinos,  S.  C.  and  W.  C.  Reynolds  (1994).  A  structure-based  model  for  the  rapid  dis¬ 
tortion  of  homogeneous  turbulence.  Technical  Report  TF-61,  Mechanical  Engineering 
Department,  Stanford  University 

Kassinos,  S.  C.,  W.  C.  Reynolds,  and  M.  M.  Rogers  (2001).  One-point  turbulence  struc¬ 
ture  tensors.  J.  Fluid  Mech.  428,  213-248. 

Kim,  J.  (1992).  Turbulence  structure  database  in  channel  flow  at  Rer  =  385.  Unpublished. 

Mansour,  N.  N.,  T.-H.  Shih,  and  W.  C.  Reynolds  (1991).  The  effects  of  rotation  on 
initially  anisotropic  homogeneous  flows.  Phys.  Fluids  5(10),  2421-2425. 

Moser,  R.  D.,  J.  Kim,  and  N.  N.  Mansour  (1999).  Direct  numerical  simulation  of  turbulent 
channel  flow  up  to  Rer  —  590.  Phys.  Fluids  11  ( 4),  943-945. 

Oberlack,  M.,  W.  Cabot,  and  M.  M.  Rogers  (1998).  Group  analysis,  DNS  and  modeling 
of  turbulent  channel  flow  with  streamwise  rotation.  In  Studying  Turbulence  Using 
Numerical  Simulation  Databases  VII.  Proceedings  of  the  1998  Summer  Program ,  pp. 
221-242.  NASA  Ames/Stanford  Center  for  Turbulence  Research. 

Obukhov,  A.  M.  (1941).  On  the  energy  distribution  in  the  spectrum  of  a  turbulent  flow. 
DokL  Akad.  Nauk  SSSR  32 ,  22-24.  As  reprinted  in  C.  R.  (Dokl.)  Acad.  Sci.  URSS 
32,  19-21  (1963). 

Osterlund,  J.  M.,  A.  V.  Johansson,  H.  M.  Nagib,  and  M.  H.  Hites  (2000).  A  note  on  the 
overlap  region  in  turbulent  boundary  layers.  Phys.  Fluids  12(1),  1-4. 

Parneix,  S.,  P.  A.  Durbin,  and  M.  Behnia  (1998).  Computation  of  3-D  turbulent  boundary 
layers  using  the  V2F  model.  Flow ,  Turbulence  and  Combustion  60,  19-46. 

Paullay,  A.  J.,  R.  E.  Melnik,  A.  Rubel,  S.  Rudman,  and  M.  J.  Siclari  (1985).  Similarity 
solutions  for  plane  and  radial  jets  using  a  h-e  turbulence  model.  J.  Fluids  Engg.  107, 
79-85. 

Pettersson  Reif,  B.  A.,  P.  A.  Durbin,  and  A.  Ooi  (1999).  Modeling  rotational  effects  in 
eddy-viscosity  closures.  Int.  J.  Heat  Fluid  Flow  20,  563-573. 

Reynolds,  W.  C.  (1991).  Towards  a  structure-based  turbulence  model.  In  T.  B.  Gatski, 
S.  Sarkar,  and  C.  G.  Speziale  (Eds.),  Studies  in  Turbulence,  Lumley  60th  birthday 
Symposium,  pp.  76-80.  Springer-Verlag,  New  York. 

Reynolds,  W.  C.  and  A.  K.  M.  F.  Hussain  (1972).  The  mechanics  of  an  organized  wave  in 


REFERENCES 


117 


turbulent  shear  flow.  Part  3.  Theoretical  models  and  comparisons  with  experiments. 
J.  Fluid  Mech.  54(2),  263-288. 

Reynolds,  W.  C.  and  S.  C.  Kassinos  (1995).  One-point  modelling  of  rapidly  deformed 
homogeneous  turbulence.  Proc.  Roy .  Soc.  London  A  451,  87-104.  Osborne  Reynolds 
Centenary  Volume. 

Reynolds,  W.  C.,  S.  C.  Kassinos,  C.  A.  Langer,  and  S.  L.  Haire  (2000).  New  directions  in 
turbulence  modeling.  Presented  at  the  Third  International  Symposium  on  Turbulence , 
Heat ,  and  Mass  Transfer.  Nagoya,  Japan  April  3-6. 

Reynolds,  W.  C.,  C.  A.  Langer,  and  S.  C.  Kassinos  (2002).  Structure  and  scales  in 
turbulence  modeling.  Phys.  Fluids  14( 7),  2485-2492. 

Ristorcelli,  J.  R.,  J.  L.  Lumley,  and  R.  Abid  (1995).  A  rapid-pressure  covariance  repre¬ 
sentation  consistent  with  the  Taylor-Proudman  theorem  materially  frame  indifferent 
in  the  two-dimensional  limit.  J.  Fluid  Mech.  292 ,  111-152. 

Robinson,  D.  F.,  J.  E.  Harris,  and  H.  A.  Hassan  (1995).  Unified  turbulence  closure  model 
for  axisymmetric  and  planar  free  shear  flows.  AIAA  J.  3S(Y2),  2325-2331. 

Singh,  K.  and  J.  L.  Lumley  (1971).  Extension  of  Heisenberg’s  model  for  spectral  transfer 
to  second-order  fluids  in  turbulent  flow.  Appl  Sci.  Res.  24  (2-3),  187-202. 

Speziale,  C.  G.  (1981).  Some  interesting  properties  of  two-dimensional  turbulence.  Phys. 
Fluids  24  ( 8),  1425-1427. 

Speziale,  C.  G.  and  N.  Mac  Giolla  Mhuiris  (1989).  On  the  prediction  of  equilibrium  states 
in  homogeneous  turbulence.  J.  Fluid  Mech.  209 ,  591-615. 

Speziale,  C.  G.,  S.  Sarkar,  and  T.  B.  Gatski  (1991).  Modelling  the  pressure-strain  corre¬ 
lation  of  turbulence  -  an  invariant  dynamical  systems  approach.  J.  Fluid  Mech.  227 , 
245-272. 

Strang,  G.  (1993).  Introduction  to  Linear  Algebra.  WeUesley-Cambridge  Press. 


